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PREFACE 


The  work  reported  here  was  performed  in  the  Optical  and  IR  Science  Laboratory  of 
the  Advanced  Concepts  Division  (which  is  now  the  Electro-Optics  Laboratory),  the 
Environmental  Research  Institute  of  Michigan  (ERIM).  The  work  was  sponsored  by  the 
Office  of  Naval  Research  (ONR),  Boston,  Contract  No.  N00014-93-C-0093,  funded 
from  the  Innovative  Science  and  Technology  Office  of  the  Strategic  Defense  Initiative 
Office  (SDIO/IST).  The  project  monitors  at  ONR  were  Dr.  Robert  Mongeon,  Mr. 
William  Miceli,  and  Dr.  William  Stachnik. 

This  final  report  covers  work  performed  from  June  1993  to  July  1994.  The 
principal  investigator  was  Dr.  James  R.  Fienup.  Additional  contributors  to  this  work 
also  included  Dr.  Jack  N.  Cederquist,  Ms.  Ann  M.  Kowalczyk,  Dr.  Joseph  C.  Marron, 
Dr.  Richard  G.  Paxman,  Mr.  Kirk  S.  Schroeder,  and  Ms.  Joy  E.  VanBuhler. 


AcoasBion  For 

EfIS  GRA&I 

DTIC  'MS 

□ 

Unannoun®  e& 

□ 

Justification — 

- - - 

By - 

Dlstgl-tntlon/ . 

Aval labil tty  Cvfijb 


veil  s»4/0F 

Sp'SC  iti.1 


®>ERIM 


(blank  page) 


IV 


®>ERIM 


TABLE  OF  CONTENTS 

PREFACE . iii 

TABLE  OF  CONTENTS . v 

LIST  OF  FIGURES . vii 

LIST  OF  TABLES . viii 

1  INTRODUCTION  AND  OVERVIEW . 1 

2  PLAN  FOR  FIELD  EXPERIMENTS . 5 

2.1  FOUR  PLANNED  EXPERIMENTS . . . 5 

2.2  SCIENCE  QUESTIONS . 8 

2.3  EXPERIMENTAL  PLAN  DETAILS . 10 

2.3.1  Speckle  Camera  Direct  Detection  Experiments . 10 

2.3.2  Short-Range  Distributed-Array  Heterodyne  Experiments . 16 

2.3.3  Three-Dimensional  Imaging  Experiments . 20 

2.3.4  Long-Range  Imaging  Experiments . . . 24 

2.4  SECTION  2  REFERENCES . 26 

3  TRADE-OFF  BETWEEN  ARRAY  SPARSITY  AND  NUMBER  OF 

REALIZATIONS . 27 

3.1  DIFFRACTION-LIMITED  COHERENT  IMAGING . 27 

3.2  IMAGING  CORRELOGRAPHY . 38 

3.3  SIMULATION  DETAILS . 41 

3.4  SECTION  3  REFERENCES . 46 

4  PHASING  OF  SPARSE  ARRAYS  OF  HETERODYNE  RECEIVERS . 47 

4.1  INTRODUCTION . 47 

4.2  IMAGING  MODEL . 47 

4.3  GENERALIZED  ITERATIVE  TRANSFORM  ALGORITHM . 48 

4.4  MAXIMIZING  SHARPNESS . 51 

4.5  COMPARISON  WITH  MAXIMIZING  SHARPNESS . 53 

4.6  IMAGE  RECONSTRUCTION  EXAMPLE . 55 

4.7  CONCLUSIONS . 56 

4.8  SECTION  4  REFERENCES . 58 

5  TURBULENCE  ISSUES  FOR  LASER-IMAGING  EXPERIMENTS . 61 

5.1  INTRODUCTION . 61 

5.2  TURBULENCE  STRENGTH . 61 

5.3  ANISOPLANATISM . 62 

5.4  RECOMMENDATIONS  AND  OUTSTANDING  ISSUES . 67 

5.5  SECTION  5  REFERENCES . 68 

6  NEW  3-D  IMAGING  MODALITY  USING  AN  OPACITY  CONSTRAINT . 69 

6.1  INTRODUCTION . 69 

6.2  IMAGING  CONCEPT . 69 

6.3  SYNTHETIC- APERTURE  DESIGN . 70 


v 


®>ERIM 


6.4  SERIAL-FREQUENCY  ILLUMINATION . 71 

6.5  PARALLEL-FREQUENCY  ILLUMINATION . 74 

6.6  CONCLUSIONS . 77 

6.7  SECTION  6  REFERENCES . . 78 

7  LABORATORY  DEMONSTRATION  OF  QUASI-HETERODYNE  IMAGING....  79 

7.1  LABORATORY  DATA  COLLECTION . 79 

7.2  IMAGE  RECONSTRUCTION . . 80 

7.3  SECTION  7  REFERENCES . . . 81 

8  DETERMINISTIC  SPARSE  ARRAY  DESIGN . 83 

8.1  SECTION  8  REFERENCES . 85 

9  INTERMEDIATE-PLANE  DETECTION . 86 

APPENDIX  A . 89 

APPENDIX  . . HI 

APPENDIX  C . H9 

APPENDIX  D . 131 


vi 


®>ERIM 


LIST  OF  FIGURES 


Figure  1-1.  Unconventional  Aperture-Plane  Imaging  with  Lasers . 1 

Figure  2-1.  Illustration  of  the  Relations  Among  the  Major  Sets  of  Experiments 

(boxes) . 5 

Figure  2-2.  Speckle  Camera  for  Direct  Detection  Imaging  Experiments . 1 1 

Figure  3-1.  Imaging  Four-Points  Object  Through  a  Sparse  Aperture  of  21  Circles...  30 

Figure  3-2.  Transfer  Functions  for  Imaging  Through  a  Sparse  Aperture  of  21 

Circles . 33 

Figure  3-3.  Imaging  Bus+Rv  Object  Through  a  Sparse  Aperture  of  21  Circles . 34 

Figure  3-4.  Spatial  Power  Spectra  of  Three  Objects . 35 

Figure  3-5.  Imaging  Through  a  22-Element  Detector  Array . 37 

Figure  3-6.  Fourier  Magnitude  and  Reconstructed  Image  from  Correlography . 39 

Figure  3-7.  Image  from  Correlography  for  Different  Numbers  of  Realizations . 40 

Figure  3-8.  Model  of  Signal-to-Noise  Ratio  Versus  Spatial  Frequency . 44 

Figure  3-9.  Model  of  Noise-to-Signal  Ratio  Versus  Spatial  Frequency . 45 

Figure  3-10.  Model  of  Signal-to-Noise  Ratio  Versus  Spatial  Frequency . 45 

Figure  3-11.  Model  of  Noise-to-Signal  Ratio  Versus  Spatial  Frequency . 46 

Figure  4-1.  Computer  Simulation  of  Sparse  Heterodyne  Array  Phasing  by  the 

Iterative  Transform  Algorithm . 57 

9 

Figure  5-1.  Cn  Measurements . 63 

Figure  5-2.  RMS  Phase  Difference  as  a  Function  of  Field  Position . 65 

Figure  6-1.  Detector  Geometry . 72 

Figure  6-2.  Time-Frequency  Representation  of  Step-Frequency  Illumination  Pulse  ..  73 

Figure  6-3.  Parallel-Frequency  Illumination  Mode . 75 

Figure  7-1.  Experimental  Setup  Used  for  Laboratory  Data  Collection . 79 

Figure  7-2.  Image  from  Quasi-Heterodyne  Laboratory  Experiment . 81 

Figure  8-1.  U  Array  Pattern . 84 

Figure  8-2.  Arrow  Array  Pattern . 84 

Figure  8-3.  Comparison  of  Number  of  Subapertures  Needed  to  Synthesize  a  Filled 

OTF  as  a  Function  of  OTF  (Autocorrelation)  Size . 85 

Figure  8-4.  Box+2  Array  Pattern . 85 

Figure  9-1.  Detector  Array  in  a  Plane  Behind  the  Focal  Plane  of  the  Lens . 87 


Vll 


LIST  OF  TABLES 


Table  2-1.  Impact  of  Various  Parameters  on  Speckle  Camera  Experiments . 13 

Table  5-1.  Correlation  Diameter,  r0,  as  a  Function  of  Cn  and  Range,  R . 62 

Table  5-2.  Isoplanatic  Patch  Diameter  as  a  Function  of  Turbulence  Strength,  C^2, 

and  Range,  R . 66 

Table  5-3.  Aperture  Diameter  as  a  Function  of  Turbulence  Strength,  Cn  . 67 

Table  8-1.  Number  of  Subapertures  to  Fill  (2L-l)x(2L-l)  Autocorrelation  Space....  84 


vm 


®>ERIM 


1  INTRODUCTION  AND  OVERVIEW 


This  report  describes  an  effort  to  examine  a  variety  of  unconventional  imaging 
experiments  that  could  be  carried  out  at  the  Innovative  Science  and  Technology 
Experimental  Facility  (ISTEF)  at  Cape  Canaveral,  FL.  These  unconventional  imaging 
modalities  all  share  the  following  attributes:  they  are  based  on  illuminating  a  target  by  a 
laser  and  detecting  the  reflected  light  in  a  non-imaged  (pupil  or  aperture)  plane,  as 
illustrated  in  Figure  1-1.  This  allows  for  a  very  light-weight,  rapidly  steerable, 
inexpensive  receiver  capable  of  very  fine  resolution.  Such  imaging  approaches  could 
have  a  large  impact  for  such  applications  as  theater  missile  defense,  interceptor  updates, 
surveillance,  and  medical  imaging.  The  imaging  approaches  we  studied  include 
heterodyne  array  imaging,  quasi-heterodyne  imaging,  laser  correlography,  and  a  new 
imaging  modality:  3-D  imaging  by  phase  retrieval  using  an  opacity  constraint. 
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o°o°o  °  Sparse  array 

Atmosphere/0  o  o  =  aff°rdable, 
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Figure  1-1.  Unconventional  Aperture-Plane  Imaging  with  Lasers. 


Under  this  effort  we  accomplished  the  following. 
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We  planned  sparse-array  unconventional  imaging  experiments  for  ISTEF.  We 
outlined  experimental  considerations  for  experiments  to  demonstrate  imaging  with 
sparse  arrays  of  heterodyne  receivers.  Our  plans  are  described  in  Section  2. 

Through  computer  simulations  we  demonstrated  the  trade-off  between  the  number 
of  sparse-array  subapertures  and  the  number  of  speckle  realizations  for  heterodyne 
array  imaging  and  imaging  correlography.  We  showed  that  when  imaging  a  low- 
contrast,  extended  target  with  sparse  arrays,  many  dozens  or  hundreds  of  speckle 
realizations  are  required  for  good-quality  imagery,  whereas  for  a  target  consisting  of 
only  a  few  bright  points,  only  a  few  speckle  realizations  are  required.  This  is  described 
in  Section  3. 

We  developed  two  new  algorithms  for  phasing  sparse  arrays  of  heterodyne 
receivers,  which  no  previous  algorithms  could  do  without  additional  information. 
These  algorithms  would  correct  for  phase  errors  due  to  both  heterodyne-receiver  local- 
oscillator  path-length  errors  and  atmospheric  turbulence.  One  of  the  algorithms  is 
based  on  the  iterative  transform  algorithm  using  a  support  constraint,  and  the  other  is 
an  image  sharpening  algorithm.  Section  4  gives  derivations  of  the  algorithms  and  shows 
successful  phase-error  correction  results  on  computer-simulated  data. 

We  performed  analysis  of  near-ground  atmospheric  turbulence  and  derived  limits  of 
field-of-view  and  detector  sizes  for  the  1  km  and  10  km  ranges  at  ISTEF.  We  found  the 
most  limiting  factor  to  be  anisoplanatism.  This  is  described  in  Section  5. 

We  analyzed  a  multiple-wavelength  data-collection  system  for  use  with  3-D  phase 
retrieval  with  an  opacity  constraint.  We  performed  computer  simulations  of  3-D 
(angle-angle-wavelength)  laser  intensity  data,  developed  a  phase  retrieval  algorithm 
using  an  opacity  constraint,  and  demonstrated  image  reconstruction  from  the  simulated 
data.  This  represents  the  first  demonstration  of  this  novel  3-D  imaging  concept.  This  is 
described  in  Section  6. 
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We  set  up  a  laboratory  experiment  and  gathered  quasi-heterodyne  data  with  multiple 
speckle  realizations  and  successfully  processed  the  data  into  an  image.  This  is  described 
in  Section  7. 

We  invented  a  new,  efficient,  deterministic  sparse-array  design,  called  the  “arrow,” 
which  fills  optical  transfer  function  (OTF)  space  with  a  very  small  number  of 
subapertures  and  is  scalable  to  any  arbitrary  array  size.  It  and  a  companion  that  gives 
double  redundancy  in  OTF  space  are  described  in  Section  8. 

We  analyzed  a  novel  sensing  method  for  gathering  data  with  a  CCD  array  that  is  in 
an  intermediate  plane  that  is  in  neither  a  focal  plane  nor  in  a  plane  conjugate  to  the 
pupil  plane.  This  allows  us  to  perform  pupil-plane  detection  while  simultaneously 
obtaining  some  background  suppression  by  limiting  the  field-of-view,  and  optimally 
trading  off  these  two  aspects. 

To  make  this  report  more  self-contained  and  to  avoid  the  need  for  extensive  review 
material  within  the  main  body  of  the  report,  we  have  included,  as  Appendix  A,  a 
description  the  theory  behind  most  of  the  imaging  modalities  analyzed  in  this  report.  It 
is  from  J.R.  Fienup,  “Unconventional  Imaging  Systems:  Image  Formation  from  Non- 
imaged  Laser  Speckle  Patterns,”  in  Emerging  Systems  and  Technologies.  S.R. 
Robinson,  Ed.,  Vol.  8  of  The  Infrared  &  Electro-optical  Systems  Handbook.  (ERIM, 
Ann  Arbor,  1993),  Chapter  1.5,  pp.92-109.  A  reader  unfamiliar  with  pupil-plane  laser 
imaging  modalities  would  benefit  from  reading  Appendix  A  first. 

References  are  at  the  end  of  each  section. 

Portions  of  this  effort  are  described  in  the  following  presentations  and  publications. 
We  have  included  as  Appendices  B,  C,  and  D  the  first  three,  which  have  appeared  as 
proceedings  papers. 

J.R.  Fienup,  A.M.  Kowalczyk  and  J.E.  Van  Buhler,  “Phasing  Sparse  Arrays 
of  Heterodyne  Receivers,”  presented  in  Orlando,  Fla.,  and  published  in  Proc. 
SPIE  2241-15.  Inverse  Optics  III  (April  1994),  pp.  127-131. 
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R.G.  Paxman,  J.H.  Seldin,  J.R.  Fienup,  and  J.C.  Marron,  “Use  of  an  Opacity 
Constraint  in  Three-Dimensional  Imaging,”  presented  in  Orlando,  FL,  and 
published  in  Proc.  SPIE  2241-14.  Inverse  Optics  in  (April  1994),  pp.  116- 
126. 

J.C.  Marron  and  R.G.  Paxman,  “Pulsed  Heterodyne-Array  Imaging,” 
presented  in  Orlando,  FL,  and  published  in  Proc.  SPIE  2241-13,  Inverse 
Optics  III  (April  1994),  pp.  111-115. 

J.R.  Fienup,  “Phasing  Distributed  Arrays  of  Heterodyne  Receivers,”  presented 
to  the  Annual  Meeting  of  the  Optical  Society  of  America,  Dallas,  TX,  October 
1994. 

J.C.  Marron,  J.R.  Fienup,  and  R.G.  Paxman,  “Unconventional  3-D  Laser 
Imaging,”  submitted  for  presentation  to  Lasers  ’94,  Quebec,  Canada, 
December  1994. 
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2  PLAN  FOR  FIELD  EXPERIMENTS 


We  envision  a  set  of  four  major  experimental  thrusts  to  be  performed  at  ISTEF:  (1) 
direct-detection  imaging,  (2)  short-range  distributed  heterodyne  array  imaging,  (3)  3-D 
imaging,  and  (4)  long-range  imaging.  Each  major  experiment  would  build  on  the 
infrastructure  and  knowledge  resulting  from  the  previous  experiments.  The 
relationships  among  these  experiments  are  illustrated  in  Figure  2-1.  They  are  described 
briefly  in  the  Section  2.1  and  in  more  detail  in  Section  2.3. 


Figure  2-1.  Illustration  of  the  Relations  Among  the  Major  Sets  of  Experiments  (boxes). 
The  ovals  show  how  infrastructure  capabilities  will  be  introduced  to  enable  each  new 
set  of  experiments  in  the  experimental  sequence. 


2.1  FOUR  PLANNED  EXPERIMENTS 


In  our  first  major  set  of  experiments,  we  would  reconstruct  high-quality  images 
from  speckle  data  collected  at  the  1  km  range  at  ISTEF.  These  speckle  data  would  be 
collected  using  a  speckle  camera  -  a  high-quality  detector  array  in  the  reimaged  pupil 
plane  of  a  telescope.  This  initial  sensor,  which  would  be  constructed  and  demonstrated 
in  Years  1  and  2  of  the  effort,  would  allow  us  to  demonstrate  image  formation  by 
modalities  such  as  imaging  correlography  and  coherent-imaging  by  phase  retrieval.  A 
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further  goal  of  these  experiments  is  to  determine  our  ability  to  measure  key  target 
parameters  such  as  size,  rotation  rate,  and  rotation  axis,  and  to  answer  a  variety  of 
other  science  questions.  The  activities  with  the  speckle  camera  would  go  on  throughout 
the  three  years  of  the  program. 

The  second  major  set  of  experiments  would  involve  the  demonstration  of  imaging 
with  a  distributed  array  of  heterodyne  receivers.  These  experiments  would  be 
performed  in  Years  1  and  2.  The  eight  heterodyne  receiver  elements  built  by  CREOL 
would  be  arranged  in  a  sparse  (quasi-non-redundant)  1-D  array.  Imaging  experiments 
would  be  performed  using  target  rotation  to  synthesize  the  second  dimension  of  the 
receiver  aperture.  Experiments  would  also  be  conducted  with  the  array  elements  in  a 
non-redundant  2-D  array.  We  would  conduct  these  experiments  at  the  1  km  and  3  km 
ranges  at  ISTEF  in  cooperation  with  other  contractors.  To  obtain  meaningful  imagery 
in  Year  2  experiments  with  non-synthetic  apertures,  a  32-element  array  of  receiver 
elements  would  be  built.  In  Year  1  of  the  program,  contractors  would  produce  a  high- 
performance,  cost-effective  design  that  would  be  the  key  to  an  affordable  array  of  this 
size.  These  experiments  would  make  use  of  knowledge  from  the  previous  stages  of  the 
program. 

A  third  set  of  experiments,  conducted  primarily  during  Years  2  and  3,  would  be 
designed  to  examine  a  technique  for  3-D  imaging.  We  refer  to  the  3-D  imaging  method 
as  Phase  Retrieval  with  an  Opacity  Constraint  for  LAser  IMaging  (PROCLAIM).  This 
technique  relies  on  measurement  of  the  speckle  patterns  recorded  at  a  series  of 
wavelengths.  A  tunable  illumination  source  would  be  added  to  the  ISTEF  experimental 
suite  to  perform  these  experiments.  The  speckle  camera  would  be  used  to  record 
speckle  patterns  at  a  series  of  wavelengths.  An  algorithm  invented  at  ERIM  that  makes 
use  of  an  opacity  constraint  would  be  tested  for  the  reconstruction  of  3-D  images. 
These  experiments  would  build  upon  the  experience  gained  in  the  speckle  camera 
experiments. 

Long-range  imaging  experiments  would  comprise  the  final  set  of  experiments  in  this 
program.  These  experiments  would  be  begun  during  Year  1  and  performed  primarily 
during  Years  2  and  3.  They  would  use  the  imaging  modalities  developed  earlier  in  the 
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program,  with  emphasis  on  the  32-element  array  of  heterodyne  receivers.  Carefully 
controlled  experiments  would  be  designed  to  perform  imaging  at  ranges  of  from  10  to 
1CI0  km.  Potential  targets  include  tethered  balloons,  manned  or  unmanned  aircraft 
flying  dedicated  paths  over  the  ISTEF  range,  and  sounding  rocket  borne  payloads,  such 
as  balloons.  Each  of  the  receivers  developed  during  the  previous  phases  of  the  program 
would  be  modified  and  used  for  these  experiments. 

At  each  stage,  the  imaging  testbed  would  build  upon  the  previous  stages,  both 
conceptually  and  using,  to  the  extent  possible,  previously  acquired  hardware. 
Furthermore,  at  each  stage  we  would  be  exploring  science  questions  (described  earlier), 
and  so  each  year's  effort  would  result  in  tangible  accomplishments  as  well  as  building 
toward  the  further  stages 

The  ultimate  goal  of  the  experiments  is  to  image  distant  targets.  However,  for  the 
less-mature  direct-detection  imaging  modalities  we  would  start  with  receiver 
characterization  locally  and  then  experiments  indoors  at  ISTEF  before  moving  to  the 
outside  range.  Then  as  the  imaging  modalities  mature,  we  would  perform  experiments 
on  the  ISTEF  outdoor  1  km  range  and  later  employ  longer  up-looking  paths  using 
manned  or  unmanned  aircraft,  weather  balloons,  or  sounding-rocket  payloads  as 
targets. 

In  parallel  with  the  experimental  effort,  we  would  perform  extensive  analysis, 
computer  simulations,  and  algorithm  development.  We  would  predict  received  light 
levels,  calculate  signal-to-noise-ratios,  and  derive  information-theoretic  Cramer-Rao 
lower  bounds  on  performance.  We  would  also  analyze  the  advantages  and  disadvantages 
of  sparse  arrays  over  filled  arrays.  Simulations  would  include  laser  power  and 
coherence  properties,  solid  modeling  of  three-dimensional  targets,  modeling  target 
surface  scattering,  polarization  effects,  atmospheric  turbulence  (including  phase  errors, 
scintillation,  and  anisoplanatism),  receiver  aperture  configuration,  background  light, 
detector  area  and  efficiency,  and  receiver  noise.  We  would  reconstruct  images  from 
simulated  data  and  measure  their  quality  (rms  error,  resolution)  to  predict  the 
performance  of  any  given  imaging  scenario.  The  analyses  and  simulations  would  define 
the  expected  performance  of  the  candidate  imaging  modalities.  The  results  of  the 
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analysis  and  simulation  efforts  would  be  readily  applicable  to  a  wide  range  of  imaging 
modalities  and  scenarios.  The  experiments  would  provide  data  points  that  would  either 
confirm  the  analyses  and  simulations  or  cause  us  to  refine  them  so  that  we  can  more 
accurately  predict  actual  performance.  Algorithm  development  would  include 
optimizing  sparse-array  design  algorithms,  phasing-up  algorithms,  and  phase  retrieval 
algorithms,  and  investigating  the  use  of  new  approaches  such  as  genetic  algorithms. 

We  also  considered  other  possible  sets  of  experiments,  including  sparse-array  direct- 
detection  experiments  and  quasi- heterodyne  experiments.  For  both  of  these  we  decided 
that,  given  budgetary  considerations,  they  had  lower  priority  and  would  yield  less 
significant  results  than  the  four  experiments  outlined  above.  Both  of  these  experiments 
would  be  beneficial  as  precursors  to  the  sparse  heterodyne  array  experiments. 
However,  since  the  sparse  heterodyne  array  experiments  could  take  place  at  an  early 
date,  there  is  less  need  for  the  sparse-array  direct-detection  experiments.  The  sparse 
array  of  heterodyne  receivers  can  provide  similar  data  as  the  sparse-array  direct- 
detection  experiments  by  simply  discarding  the  phase  information.  Similarly,  the 
speckle-camera  experiments  can  provide  similar  data  as  the  sparse-array  direct- 
detection  experiments  by  simply  discarding  the  information  from  most  of  pixels  in  the 
CCD  array.  In  that  way  the  speckle  camera  can  explore  the  effect  of  many  different 
sparse  array  patterns  from  a  single  filled  array  of  data.  Many  aspects  of  the  quasi¬ 
heterodyne  array  experiments  can  also  be  covered  by  either  the  speckle  camera 
experiments  or  the  sparse  heterodyne  array  experiments.  However,  the  quasi¬ 
heterodyne  experiments  would  have  allowed  for  testing  a  variety  of  sparse  array 
patterns  in  a  much  more  efficient  manner.  From  one  set  of  filled-aperture  quasi¬ 
heterodyne  data  we  could,  by  processing  selected  groups  of  pixels  from  the  array,  try 
out  many  different  sparse-array  patterns. 

Section  2.3  give  additional  details  of  the  experimental  plan. 

2.2  SCIENCE  QUESTIONS 

Much  of  the  effort  in  the  planned  experiments  would  entail  the  fabrication, 
assembly,  and  testing  of  the  optical,  mechanical,  and  electronic  equipment,  integrating 
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it  at  the  ISTEF  site,  and  performing  imaging  experiments.  In  addition,  we  would 
conduct  experiments  that  would  address  numerous  scientific  questions.  Some  of  them 
are  as  follows. 

Imaging  Performance  The  most  important  question  is:  what  is  the  performance  of  a 
variety  of  active  imaging  modalities  employing  sparse,  distributed  arrays  of  detectors. 
We  would  determine  the  trade-off  between  array  sparsity,  number  of  speckle 
realizations,  and  laser  power. 

Target  Effects  What  are  the  effects  of  the  target  on  the  received  speckle  pattern  and  on 
the  image?  In  particular,  what  is  the  effect  of  target  flexure  and  of  multiple  targets  in 
the  field  of  view  rotating  at  different  angular  rates.  What  is  the  effect  of  target  size, 
shape,  and  texture  on  its  spatial  power  spectrum  and  on  the  quality  of  an  image  of  it? 
What  is  the  effect  of  the  target  geometry  and  surface  microstructure  on  polarization? 

Target  Parameter  Estimation  Short  of  reconstructing  an  image,  how  well  can  we 
determine  target  diameter  and  other  target  characteristics  from  the  aperture-plane 
speckles?  How  can  we  optimally  estimate  target  motion  from  the  motion  of  the  detected 
speckles? 

Atmospheric  Effects  How  does  imaging  performance  correlate  with  the  degree  of 
atmospheric  turbulence  (for  both  the  illumination  and  the  backscattered  light)?  What  is 
the  effect  of  anisoplanatism  on  imagery,  and  how  can  it  be  overcome?  How  much  must 
we  reduce  the  back-scattered  light  from  aerosol  scattering,  and  to  what  degree  do 
various  approaches  (such  as  using  a  finite  bistatic  angle,  range  gating,  spectral  filtering 
and  spatial  filtering)  or  a  combination  of  approaches  reduce  it?  What  laser  pulse 
duration  is  needed  to  adequately  freeze  the  atmosphere  for  different  imaging  paths? 

Illuminator  What  is  the  effect  of  partial  coherence  of  the  laser?  What  wave  forms  are 
optimal  for  the  different  imaging  modalities? 

Array  Pattern  What  is  the  effect  of  the  configuration  of  a  sparse  array?  In  particular, 
how  sparse  can  the  receiver  array  be  while  still  achieving  good  image  quality?  How 
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important  is  it  that  the  autocorrelation  of  the  array  pattern  fill  Fourier  (“u-v”)  space? 
What  spatial  frequencies  are  most  critical  to  preserve?  What  geometric  properties  are 
desirable  in  the  configuration  of  a  sparse  array? 

Other  For  operation  during  daylight  hours,  how  much  suppression  of  sunlight  reflected 
from  the  object,  and  scattered  in  the  atmosphere,  is  needed?  What  is  the  value  of  a  3-D 
signature  relative  to  a  2-D  signature? 

2.3  EXPERIMENTAL  PLAN  DETAILS 

In  what  follows  we  give  more  detailed  plans  for  the  experiments  at  ISTEF  that  were 
outlined  in  Section  2.1. 

2.3.1  Speckle  Camera  Direct  Detection  Experiments 

We  suggest  a  program  that  would  begin  with  direct  detection  experiments.  Direct 
detection  in  the  aperture  plane  has  two  primary  features  that  make  it  attractive  for 
imaging  experiments.  First,  it  is  several  times  less  expensive  to  perform  direct 
detection  imaging  than  to  use  an  array  of  heterodyne  receivers  and,  second,  direct 
detection  imaging  is  not  as  sensitive  to  atmospheric  phase  errors  because  intensity 
detection  is  employed.  In  addition,  direct  detection  offers  flexibility  for  proving 
imaging  concepts,  evaluating  performance  versus  parameters  such  as  aperture  fill,  and 
analyzing  science  issues.  To  support  the  development  of  unconventional  imaging,  we 
would  construct  a  speckle  camera  and  conduct  direct  detection  experiments  with  it.  In 
the  sections  below  we  give  a  program  approach  for  this  activity. 

2.3. 1.1  Speckle  Camera  Experiments 

The  primary  motivation  for  the  speckle  camera  experiments  is  to  produce  a  low- 
risk,  low-cost  sensor  for  early  demonstration  of  unconventional  imaging  on  the  outdoor 
test  range  at  ISTEF.  To  accomplish  this,  we  would  construct  a  speckle  camera 
experiment  as  shown  in  Figure  2-2.  The  components  of  this  experiment  include  a  CW 
laser  for  illumination  of  targets,  target  rotation  stage  for  generating  multiple 
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realizations  of  speckle  patterns,  a  high-speed  CCD  detector  array,  frame  grabber, 
computer  and  transmit/receive  optics.  This  camera  would  be  able  to  sample  speckle 
patterns  from  the  object  at  frame  rates  up  to  the  KHz  range.  Our  approach  is  to  use  a 
class  of  low-cost,  digital-output  CCD  cameras  manufactured  by  Dalsa  Inc.  These 
cameras  are  available  in  several  array  formats  (32  x  32,  64  x  64,  128  x  128,  and 
256  x  256)  and  operate  at  maximum  pixel  rates  of  16  MHz.  The  output  of  the  speckle 
camera  is  captured  by  a  digital-input  frame-grabber  and  then  stored  by  a  computer. 
Processing  of  the  speckle  data  can  be  performed  on  either  a  local  computer  for  initial 
processing  or  on  other  computers  to  which  the  data  can  be  transferred.  The  equipment 
would  be  tested  first  on  an  indoor  test  range  at  ERIM  and  then  on  an  outdoor  test  range 
at  ISTEF  where  long  horizontal  paths  are  available. 


Figure  2-2.  Speckle  Camera  for  Direct  Detection  Imaging  Experiments. 


The  main  experiments  to  be  conducted  using  the  speckle  camera  include  image 
formation  experiments  and  supporting  experiments  to  address  various  science 
questions.  The  two  primary  image  formation  methods  that  would  be  used  for  image 
recovery  are  imaging  correlography  [2.1],  which  yields  an  incoherent  (non- speckled) 
image,  and  complex  phase  retrieval  [2.2]  which  yields  a  coherent  (speckled  and 
complex-valued)  image.  These  methods  for  image  recovery  are  discussed  in  greater 
depth  in  Appendix  A.  The  collected  data  would  consist  of  a  filled  array  of 
measurements;  we  can  also  simulate  the  effect  of  a  sparse,  distributed  array  by 
processing  only  a  selected  subset  of  the  data.  The  science  issues  to  be  addressed  are 
discussed  in  Section  2.2.  They  include  image  quality,  the  effects  of  target  properties 
and  motion  on  the  speckle  pattern,  polarization  effects,  array-pattern  trade-offs. 
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atmospheric  effects,  and  so  forth.  The  path  at  the  outdoor  test  range  at  ISTEF  can 
provide  challenging  atmospheric  conditions  due  to  its  proximity  to  the  ground  and 
extended  path  lengths.  While  direct  detection  is  much  less  sensitive  to  turbulence  than 
heterodyne  detection,  direct  detection  is  sensitive  to  anisoplanatism,  beam  wander,  and 
break-up  (scintillation)  for  both  the  outward  and  return  paths.  As  part  of  these 
experiments,  atmospheric  conditions  would  be  monitored,  and  modifications  to  the 
beam  path,  such  as  elevating  it  from  the  ground  and  adapting  the  length,  would  be 
considered. 


Table  2-1  contains  a  discussion  of  various  parameters  and  how  they  affect  the 
speckle  camera  experiments.  As  an  example  we  give  parameters  for  a  relatively  benign 
value  of  the  atmospheric  structure  function,  Cn2.  More  stressing  turbulence  would 
cause  us  to  use  smaller  targets  and  larger  speckle  demagnifications.  Measurements  of 
Cn2  at  the  site  would  be  necessary  for  the  final  design  of  the  experiment. 
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Table  2-1.  Impact  of  Various  Parameters  on  Speckle  Camera  Experiments 


Experimental 

Parameter 

Comment 

Target  Size 

Object  size  must  be  smaller  than  isoplanatic  patch  size. 
Preliminary  analysis  gives  object  size  of  10  cm  at  1  km  range. 

Target  Material 

Variety  of  materials  would  be  evaluated.  Examples  include: 
diffuse  paints,  depolarizing  coatings,  specular  metal,  or  retro- 
paint  for  boosting  signal. 

Speckle  Size 

For  target  diameter  D  =  10  cm  at  range  R  =  1  km  and 

X  =  0.5  |im,  speckle  size  (XR/D)  is  5  mm. 

Target  Rotation 

Target  rotation  causes  speckle  patterns  to  translate.  Generation 
of  independent  speckle  realizations  requires  rotation  of  NX/D 
radians,  where  N  is  #  of  speckles  across  array.  For  N  =  32, 
rotation  of  0.16  mrad  is  required.  Measurements  of  rotation 
rate  can  be  made  for  rates  up  to  roughly  1  rad/sec. 

Path  Geometry 

Bistatic  path  is  desirable  for  background  rejection. 

Background  rejection  can  also  be  achieved  by  locating  laser 
near  target  and/or  folding  path.  Path  can  also  be  elevated  to 
reduce  turbulence. 

Signal  Level 

Signal  levels  for  diffuse  targets  are  marginal  for  highest  frame 
rates.  To  boost  signal  level,  retro-paint  can  be  used,  or 
detector-array  master  clock  rate  can  be  reduced. 

Daytime  Operation 

MODTRAN  calculations  show  that  daytime  operation  is 
feasible  with  10  nm  spectral  filter. 

Optics 

Size/Magnification 

Size  of  receiver  fore-optic  encompasses  N  x  N  speckles.  For 

N  =  32,  fore-optic  is  0.16  m  in  diameter.  For  2x2 
sampling  of  speckles  onto  16  [I  pixels,  demagnification  of 
pupil  by  156x  is  required.  Can  be  accomplished  using  a 
telescope/microscope  objective  pair. 

2.3.1.2  Algorithm  Development 


The  main  aspect  of  algorithm  development  would  be  the  application  of  previously 
demonstrated  image  recovery  algorithms  to  the  data  generated  by  the  speckle  camera. 
For  complex  phase  retrieval  it  would  also  be  necessary  to  develop  improved  estimates 
of  the  support  on  the  autocorrelation  from  the  computed  autocorrelation  function  [2.3]. 
This  is  made  difficult  because  the  autocorrelation  function  is  speckled.  For  the  case  of 
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sparse  arrays,  the  problem  is  compounded.  Algorithms  would  be  developed  to  format 
and  calibrate  data,  then  adapt  image  formation  and  speckle  analysis  software  to  the 
available  computer  platforms.  The  speckle  analysis  software  includes  features  such  as 
estimation  of  target  rotation  rates  and  rotation  axis  from  translating  speckles.  In 
addition,  control  and  diagnostic  software  would  be  required  for  (a)  capturing  multiple 
realizations  of  the  speckle  data  as  the  target  is  rotated  under  computer  control,  and  (b) 
assessing  the  quality  of  the  speckle  patterns  for  various  atmospheric  conditions  and 
detection  circumstances,  such  as  noise  level  and  sampling  geometry. 

2.3. 1.3  Experiment  Simulation 

Assessment  of  the  operation  of  image  recovery  algorithms  developed  for  the  speckle 
camera  requires  the  ability  to  simulate  speckle  patterns  that  correspond  to  the  type  of 
targets  being  imaged  under  appropriate  atmospheric  and  detection  conditions.  To 
accomplish  this,  we  would  generate  reflectance  models  of  the  targets  and  then  apply 
statistical  processes  associated  with  scattering  and  propagation  through  the  atmosphere. 
We  would  then  apply  detector  sampling  and  noise  processes  to  generate  versions  of  the 
detected  speckle  patterns.  Resulting  speckle  patterns  would  then  be  used  as  input  to 
image  recovery  and  speckle  analysis  algorithms.  Simulated  recovered  images  and 
speckle  parameters  would  then  be  compared  to  actual  images  and  parameters 
reconstructed  from  the  real  data.  This  comparison  is  valuable  for  analyzing 
nonintuitive  effects  that  may  arise  in  either  collected  or  simulated  data. 

2.3.1 .4  Experiment  Design  and  Analysis 

Early  in  the  program  we  would  conduct  an  experiment  design  and  analysis  task.  We 
have  already  performed  feasibility  calculations  of  atmospheric  turbulence  parameters 
and  signal  power  levels.  This  analysis  would  be  refined  and  used  for  detailed 
experimental  design,  which  would  allow  us  to  determine  experimental  parameters  such 
as  maximum  target  size  and  maximum  range,  and  specify  optical  components. 
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2.3. 1.5  Hardware  Fabrication/Acquisition 


Upon  completion  of  the  design  task,  the  speckle  camera  would  be  fabricated.  This 
involves  the  acquisition  of  new  hardware  components  and  fabrication  of  some 
components  such  as  mounting  hardware  for  optics. 

2.3. 1.6  Hardware  Integration/Lab  Testing 

Following  completion  of  the  hardware  fabrication  and  acquisition,  we  would 
integrate  and  test  the  speckle  camera  at  a  20-meter  indoor  test  range  at  ERIM  in  Ann 
Arbor.  At  ERIM  we  have  access  to  auxiliary  equipment,  including  electronic  test 
equipment,  that  would  be  required  for  integration  and  testing.  For  testing  at  ERIM  we 
would  use  an  existing  argon  laser  for  illuminating  targets.  This  testing  phase  would 
include  characterization  of  the  recorded  data  by  evaluating  properties  such  as  the 
speckle  statistics  and  demonstration  of  image  recovery. 


2.3. 1.7  Field  Experiment  Execution 


Once  lab  testing  is  complete,  we  would  transport  the  speckle  camera  to  the  ISTEF 
facility  for  experiment  execution.  The  laser  power  required  is  1  or  2  watts,  preferably 
in  the  green.  As  mentioned  above,  we  would  use  an  existing  argon  laser  as  a  source  at 
ERIM.  For  this  effort  we  require  an  appropriate  laser  at  ISTEF,  which  could  be  either 
an  argon  laser  or  a  frequency-doubled,  CW,  YAG  laser.  If  no  such  lasers  are  available, 
an  ERIM  laser  could  be  used.  We  would  perform  initial  experiments  at  a  relatively 
modest  range,  and  then  perform  the  main  experiments  at  the  1-km  range  and  longer 
ranges.  An  important  consideration  in  the  design  of  the  experiments  is  the  laser  power 
level.  The  desirable  received  light  from  the  target  is  significantly  above  the  noise  level 
of  the  detector.  The  energy  received  by  each  pixel  of  the  array  is  given  by 
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assuming  that  the  target  is  perfectly  diffuse  and  the  detector  size  is  equal  to  the  speckle 
size.  With  the  laser  power  PL  =  1.0  watt,  X  =  0.5  p.m,  At  =  1/12,000  second  (which 
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is  the  fastest  rate  for  the  32  x  32  CCD  camera),  object  diameter  D  =  5  cm,  object 
illumination  efficiency  T]obj=0.2,  and  path-length  transmittance  of  tit  =  0.5,  it 
follows  that  E  =  5.3  x  10'16  joules.  This  value  should  be  compared  to  the  CCD  noise 
value  of  4.48  x  10'16  joules/pixel.  Note  that  for  this  worst  case  of  short  integrations 
and  diffuse  target  the  return  is  roughly  equal  to  the  CCD  noise  level.  This  relationship 
is  improved  by  orders  of  magnitude  by  slowing  the  CCD  frame  rate  to  100  frames  per 
second  (which  gives  two  orders  of  magnitude)  and  using,  if  necessary,  retro-reflective 
paint  on  the  targets.  Also  note  that,  as  long  as  the  detector  size  scales  with  the  speckle 
size,  the  receiver  power  is  independent  of  target  range. 

2.3. 1.8  Image  Formation/Data  Analysis 

Following  speckle  camera  installation  and  data  collection,  image-formation  and 
speckle-analysis  algorithms  would  be  exercised. 

2.3.2  Short-Range  Distributed-Array  Heterodyne  Experiments 

A  significant  goal  of  the  suggested  experiments  is  to  demonstrate  imaging  with  a 
distributed  array  of  heterodyne  receivers.  To  accomplish  this,  we  would  make  use  of 
infrastructure  capabilities  developed  under  the  speckle  camera  task  described  in 
Section  2.3.1.  In  cooperation  with  other  contractors,  we  would  use  heterodyne  receiver 
elements  built  by  CREOL  to  generate  images  of  coherently  illuminated  targets. 
Challenging  aspects  of  this  task  include  designing  sparse  arrays  that  use  aperture 
synthesis  methods  to  generate  discernible  images  with  limited  numbers  of  receiver 
elements  and  applying  phase-error  correction  algorithms  that  are  required  to  phase-up 
the  receiver  array  and  to  correct  for  atmospheric  turbulence.  ERIM’s  approach  for 
these  distributed-array  experiments  makes  use  of  significant  expertise  that  we  have 
gained  through  several  years  of  research  on  unconventional  imaging  and  is  outlined 
below. 
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2.3.2. 1  Experiment  Specification 


The  distributed-array  heterodyne  experiments  have  three  primary  tasks.  First,  with 
other  contractors,  we  would  contribute  to  the  design  of  the  heterodyne  receiver 
elements  and  array  configurations.  Our  primary  input  would  be  to  establish 
performance  specifications  required  for  robust  image  formation  and  to  aid  in  the 
evaluation  of  designs  to  ensure  that  the  performance  specifications  are  met.  For 
example,  image  formation  with  heterodyne  receivers  requires  measurement  of  the 
complex-valued  optical  field  and  these  measurements  are  subject  to  errors  from  sources 
such  as  detector  noise  and  background  light.  To  correct  for  array-phasing  errors  and 
atmospheric  turbulence,  the  field  measurements  must  be  made  with  sufficient  accuracy. 
Also,  in  order  for  the  experiments  to  succeed,  certain  detector  noise  and  signal  level 
specifications  must  be  met.  Other  aspects  of  receiver  element  design  that  are  crucial  for 
experiment  success  include  the  mechanical  robustness  of  the  receiver  elements.  The 
bandwidth  of  phase  errors  imparted  by  mechanical  disturbances  of  the  optical  fiber 
network  must  be  lower  than  the  image  formation  bandwidth  to  ensure  that  the 
systematic  phase  errors  are  constant  over  the  data  collection  time.  To  accomplish  this 
task  ERIM  would  establish  performance  specifications  and  evaluate  design  approaches 
for  meeting  them  in  cooperation  with  other  contractors.  In  addition  we  would  develop 
diagnostic  software  for  the  assessment  of  measurement  quality. 

The  second  task  is  to  use  a  limited  number  of  receiver  elements  (eight)  along  with 
aperture  synthesis  methods,  to  generate  target  images  at  the  1-  and  3-km  test  ranges  at 
ISTEF.  These  experiments  would  be  conducted  in  cooperation  with  other  contractors. 
ERIM  would  contribute  to  the  areas  of  algorithm  development,  experimental  design, 
simulation,  and  image  formation.  More  specifically,  we  would  develop  experimental 
procedures  for  conducting  aperture  synthesis  based  on  object  rotation.  In  addition, 
ERIM  would  provide  algorithms  for  data  assessment,  phase  error  correction,  and 
image  formation.  The  phase-error  correction  algorithms  must  remove  phase  errors 
caused  by  array  misalignment  and  atmospheric  turbulence.  Because  the  atmospheric 
component  of  the  phase  errors  is  time-variant  it  is  likely  to  be  the  most  difficult  to 
correct,  particularly  since  aperture  synthesis  requires  an  extended  time  interval  for  data 
collection.  Ideally,  the  atmospheric  phase  would  be  frozen  (quasi- static)  during  the 
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aperture  synthesis  interval  so  the  resulting  phase  error  is  one-dimensional;  in  this  case 
the  Digital  Shearing  Laser  Interferometry  phase-error  correction  algorithm  can  be 
used  [2.4].  We  would  also  explore  ways  to  correct  phase  errors  when  they  are  changing 
within  the  data-collection  time. 

For  the  third  task  we  would  use  an  increased  number  of  receiver  elements  (32)  for 
image  generation  without  aperture  synthesis.  Key  elements  of  this  task  would  include 
the  development  of  algorithms  for  phase-error  correction  and,  in  cooperation  with 
another  contractor,  the  design  of  array  formats  that  fill  the  largest  area  with  a  small 
number  of  elements.  In  previous  work  we  have  developed  phase-up  algorithms  based 
on  image  sharpness  that  allow  for  correction  of  a  constant  phase  error  on  each  receiver 
element  that  is  attributable  to  the  array  not  being  phase  aligned  [2.5].  These  algorithms 
ideally  require  multiple  realizations  of  the  pupil-plane  speckle  field;  however,  we  have 
also  used  sharpness  to  correct  for  similar  phase  errors  encountered  in  Synthetic 
Aperture  Radar  (SAR)  when  only  a  single  realization  of  the  speckle  pattern  is  available 
[2.6].  After  phase  correction,  images  are  recovered  by  Fourier  transformation.  Image 
quality  would  then  be  assessed  by  comparison  of  reconstructed  images  to  simulated 
images. 

2A2.2  Algorithm  Development 

As  discussed  above,  a  key  element  of  these  experiments  is  algorithm  development. 
The  algorithms  of  primary  importance  apply  to  data  diagnostics,  phase  error  correction 
and  image  formation.  Once  these  algorithms  are  developed  they  would  be  ported  to  the 
primary  data  collection  and  image  formation  computers  at  ISTEF  for  application.  For 
algorithm  development  prior  to  the  availability  of  field  data,  we  require  simulated  data 
as  discussed  below. 

2.3.23  Experiment  Simulation 

Simulated  data  sets  are  required  for  algorithm  development  and  evaluation  of 
experimental  results  by  comparing  them  to  ideal  simulated  images.  ERIM  has  had 
several  years  of  experience  in  generating  simulated  imagery  including  ideal  target 
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models,  the  effects  of  atmospheric  turbulence,  imperfect  data-collection  hardware,  and 
phenomenology  modeling.  Important  aspects  of  the  simulations  required  for  this  task 
include  accurately  modeling  the  various  image-degrading  processes  (atmospheric 
turbulence,  detector  noise,  etc.)  to  evaluate  the  entire  image-formation  process  and 
diagnose  unanticipated  problems. 

2.3.2.4  Experiment  Design  and  Analysis 

Our  primary  tasks  are  to  ensure  that  the  hardware  meets  the  required  performance 
levels  for  image  formation  and  to  help  with  the  experimental  design  so  that  meaningful 
images  are  obtained  with  the  limited  number  of  receiver  elements.  To  provide  this 
support  we  would  conduct  analyses  of  the  image  formation  processes  to  determine 
various  hardware  performance  specifications  and  to  design  experiments.  In  addition, 
we  would  develop  diagnostic  software  to  use  in  the  laboratory  and  field  for  assessment 
of  data  quality. 

2.3. 2. 5  Field  Experiments 

Field  experiments  would  be  conducted  on  the  test  ranges  at  ISTEF.  During  these 
experiments  we  would  provide  support  primarily  to  ensure  that  the  data  are  of 
sufficient  quality  for  image  formation.  For  the  first  experiments,  which  use  a  small 
number  (eight)  of  receivers,  aperture  synthesis  is  recommended.  Experiments  would 
require  careful  adjustment  of  object  rotation  rate  and  careful  monitoring  of 
atmospheric  and  alignment  phase  errors.  Evaluation  of  these  effects  would  require  the 
use  of  the  diagnostic  software  for  data  assessment.  It  is  expected  that  this  data 
assessment  would  occur  as  data  are  being  collected  so  that  adjustments  can  be  promptly 
made.  The  second  set  of  experiments  would  again  require  the  careful  assessment  of  the 
data  collection  process  and  evaluation  of  output  data  using  diagnostic  software. 

2.3.2.6  Image  Formation/Data  Analysis 

We  would  exercise  our  phase  error  correction  algorithms  and  image  formation 
software,  and  to  assess  the  results  relative  to  ideal  simulated  data.  This  work  would 


19 


®>ERIM 


involve  the  two  main  experiment  types  (small  number  of  receivers  with  aperture 
synthesis  and  larger  number  of  receivers  without  aperture  synthesis).  It  is  also 
anticipated  that  several  trials  of  each  experiment  type  would  be  conducted  to  obtain 
improved  performance,  gain  experience,  and  answer  science  questions. 

2.3.3  Three-Dimensional  Imaging  Experiments 

The  PROCLAIM  imaging  modality  is  attractive  because  of  the  improved 
identification  capability  offered  by  three-dimensional  imaging  and  because  only  the 
intensity  data  are  directly  detected.  The  advantages  of  direct  detection,  including 
robustness  to  atmospheric  phase  aberrations  and  savings  in  expense  and  complexity 
relative  to  heterodyne  detection,  apply  to  the  PROCLAIM  modality.  Furthermore,  the 
speckle  camera  that  would  be  developed  under  the  first  set  of  experiments  can  be  used 
for  PROCLAIM  experiments.  The  philosophy  of  this  modality  is  to  trade  hardware 
complexity  for  increased  computational  load.  In  this  section  we  describe  the  algorithm 
development,  experiment  design,  simulations,  hardware  acquisition  and  integration,  and 
testing,  all  leading  to  a  field  experiment  in  which  the  novel  PROCLAIM  three- 
dimensional  imaging  modality  would  be  demonstrated  on  a  1-km  horizontal  test  range. 

2.3.3. 1  Three-Dimensional  Experiment  Specification 

There  are  many  illumination  strategies  that  can  be  considered  for  the  PROCLAIM 
imaging  modality.  As  discussed  in  Section  6,  we  performed  analyses  to  show  the 
preliminary  feasibility  of  using  a  serial-illumination  strategy  to  image  a  missile  or  some 
other  long-range  target  from  the  ground  with  the  PROCLAIM  concept.  For 
PROCLAIM  to  work,  the  target  should  be  relatively  stable  over  the  illumination 
sequence.  Using  the  rather  restrictive  criterion  that  the  speckle  pattern  does  not  move 
more  than  half  a  speckle  width  during  the  entire  collection  sequence,  we  arrive  at  the 
following  condition  on  the  target's  angular  frequency  of  rotation, 


where  co  is  the  object  angular  frequency,  X  is  the  wavelength,  c  is  the  speed  of  light, 
N  is  the  number  of  frequency  samples  or  equivalently  the  number  of  range  bins 
desired,  Dobj  is  the  cross-range  diameter  of  the  object,  and  Zobj  is  the  depth  of  the 
object  in  range.  Let  Dobj-  =  5,  Zobj  =  2  m,  X  =  1  |im,  and  suppose  that  we  desire 
N  -  10  range  bins.  The  maximum  angular  velocity  that  can  be  tolerated  for  such  an 
object  is  found  to  be  co  <  0.375  radians/second.  This  sample  calculation  suggests  that 
the  requirement  on  object  angular  frequency  is  not  particularly  restrictive  and  that 
there  are  likely  to  be  many  objects  that  can  be  imaged  with  serial-frequency 
illumination.  Therefore  our  ground  experiments  would  be  designed  to  investigate  this 
illumination  modality. 

The  actual  layout  of  the  experiment  would  be  closely  related  to  the  speckle  camera 
experiments,  discussed  in  Section  2.3.1.  Consequently,  the  PROCLAIM  experiments 
would  greatly  benefit  from  the  experience  gained  and  infrastructure  developed  in  the 
speckle-camera  experiments.  We  have  performed  a  preliminary  experiment  design  that 
suggests  that  a  three-dimensional  object,  such  as  a  cone  painted  with  retro-reflective 
paint,  can  be  imaged  at  a  range  of  1  km.  Anisoplanatism  is  much  more  challenging  in 
this  horizontal-imaging  scenario  than  in  the  envisioned  application  of  PROCLAIM, 
namely  looking  up  to  identify  missiles.  As  with  the  speckle  camera  experiments,  the 
object  would  need  to  be  no  larger  than  the  isoplanatic  angle,  on  the  order  of  10  cm  in 
clear  conditions.  The  same  re-imaging  optics  that  would  be  used  in  the  speckle  camera 
experiments  can  be  used  here  as  well. 

The  PROCLAIM  experiments  differ  from  the  speckle-camera  experiments  in  one 
fundamental  way:  a  frequency-tunable  laser  is  needed  for  illumination.  A  coordinated 
illumination,  frame-detection,  and  data-storage  sequence  would  be  controlled  by  a 
computer.  In  this  way  the  required  data  set,  consisting  of  far-field  speckle  (intensity) 
patterns  for  each  of  several  illuminating  frequencies,  would  be  obtained. 

2.3. 3.2  Algorithm  Development 

Since  the  PROCLAIM  concept  trades  hardware  complexity  for  increased 
computational  load,  one  might  properly  infer  that  the  algorithms  needed  for  three- 
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dimensional  image  reconstruction  are  somewhat  complicated.  A  gradient- search 
algorithm  has  been  developed  by  ERIM  researchers  and  it  has  been  demonstrated  in  a 
simple  simulation  [2.7].  This  algorithm  would  have  to  be  developed  to  accommodate 
noise  and  the  effects  of  diffraction. 

We  believe  that  the  opacity  constraint  can  be  used  not  only  for  phase  retrieval,  but 
also  for  estimating  very  tight  bounds  on  the  object  support,  given  the  object's  three- 
dimensional  autocorrelation  function  (or  equivalently,  the  object's  three-dimensional 
Fourier  intensity).  Such  a  bound  on  the  support  or  extent  of  the  object  is  called  a 
locator  set  and  ERIM  researchers  have  previously  developed  locator  set  algorithms  for 
two-dimensional  data  sets  that  can  readily  be  generalized  to  accommodate  three- 
dimensional  data  sets  [2.3].  Preliminary  analysis  suggests  that  locator  sets  for  opaque 
objects  can  be  made  very  tight  and  would  provide  excellent  initial  estimates  in  a  phase- 
retrieval  algorithm.  It  may  actually  be  possible  to  identify  a  specific  target  from  a  tight 
locator  set  that  would  give  strong  target  profile  information.  We  would  develop  both 
phase-retrieval  and  locator-set  algorithms  for  coherent  three-dimensional  opaque 
objects. 

2.3.3.3  Detailed  Experiment  Design 

This  design  would  follow  from  the  detailed  design  performed  for  the  speckle 
camera  experiments.  The  object  would  be  designed  so  that  significant  depth  information 
can  be  retrieved  and  so  that  it  can  be  rigidly  mounted.  The  detailed  design  of  the 
illumination,  detection,  and  data-storage  sequence  would  have  to  be  performed 
separately  for  the  PROCLAIM  experiments,  as  it  involves  multiple  wavelengths. 

2.3.3.4  Experiment  Simulation 

Once  the  algorithms  have  been  developed  and  the  detailed  experiment  design  has 
been  performed,  we  would  perform  simulation  experiments  that  conform  to  the 
intended  field  experiments.  The  same  speckle  simulation  algorithms  developed  in  the 
speckle-camera  experiments  would  be  used  in  this  simulation  but  with  wavelength 
diversity  added.  The  simulation  would  provide  confidence  that  we  can  succeed  with  the 
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field-experiment  demonstration.  In  addition,  the  exercise  would  provide  a  means  of 
isolating  the  source  of  non-intuitive  effects  that  may  arise  in  the  collection  of  data  or 
the  reconstruction  of  three-dimensional  objects. 


Much  of  the  hardware  required  for  the  PROCLAIM  experiments  would  already 
have  been  acquired  for  the  speckle-camera  experiments.  Of  course  the  frequency- 
tunable  laser  would  have  to  be  acquired.  ERIM  has  a  frequency-tunable  laser  that  can 
readily  be  used  in  the  laboratory  integration  and  testing  phases.  We  would  rely  on 
government-furnished  equipment  (GFE)  for  this  component  of  the  outdoor  experiment. 
The  laser  would  need  to  have  a  gain  curve  that  would  provide  the  desired  range 
resolution,  which  is  given  by 


(2.3) 


where  Av  is  the  bandwidth  over  which  the  source  is  tuned.  To  get  a  range  resolution  of 
1  cm,  we  require  Av  =  1.5  x  1010  Hz  or  A^,  =  0.25  Angstroms.  Therefore,  it  should 
be  easy  to  achieve  the  desired  range  resolution.  The  laser  would  be  tuned  with  a 
precisely  rotatable  etalon  that  is  placed  in  the  laser  cavity  and  is  electronically 
controlled  with  the  data-acquisition  computer. 


2.3. 3.6  Hardware  Integration  and  Laboratory  Testing 


Before  the  equipment  is  transported  to  the  ISTEF  test  range,  we  would  integrate  the 
experimental  hardware  in  a  laboratory  test  range  so  that  the  integrated  imaging  system 
can  be  tested  and  validated.  Again,  we  would  benefit  from  laboratory  tests  performed 
in  the  speckle-camera  experiments.  This  integration  task  would  give  us  additional 
experience  and  may  suggest  adjustments  to  the  detailed  experiment  design. 
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2.3.3.7  Execution  of  Field  Experiment 

When  the  laboratory  equipment  testing  and  integration  are  completed,  the  equipment 
would  be  transported  to  the  ISTEF  facility  to  conduct  the  field  experiment.  The 
ultimate  goal  of  the  entire  three-dimensional  imaging  task  is  to  demonstrate  the 
feasibility  and  utility  of  the  PROCLAIM  imaging  modality  on  the  ISTEF  horizontal  test 
range  at  a  range  of  1  km.  This  would  help  to  establish  the  importance  of  the  third 
dimension  of  information  for  target  discrimination.  This  could  motivate  the  extension 
of  the  heterodyne  array  imaging  from  two  dimensions  to  three. 

233Jj  Image  Formation  and  Data  Analysis 

The  phase-retrieval  algorithms  would  be  installed  on  ISTEF  computers  and  images 
would  be  formed  after  the  data  are  collected.  These  images  would  be  scored  primarily 
with  regard  to  how  well  the  target  profile  is  recovered.  Experimental  results  would  be 
compared  with  simulation  results. 

2.3.4  Long-Range  Imaging  Experiments 

The  final  set  of  experiments  to  be  conducted  in  the  series  concern  imaging  of  distant 
targets.  Emphasis  would  be  on  using  the  32-element  array  of  heterodyne  receivers 
discussed  in  Section  2.3.2,  but  experiments  with  the  speckle  camera  would  be 
conducted  as  well.  These  experiments  are  to  be  conducted  at  ranges  greater  than  3  km 
and  would  make  extensive  use  of  hardware  and  software  developed  under  previous 
tasks.  Our  approach  for  these  experiments  is  discussed  in  greater  depth  below. 

2.3.4. 1  Experiment  Specification 

The  purpose  of  these  experiments  is  to  demonstrate  long-range  imaging  using 
primarily  an  array  of  heterodyne  receivers.  However,  other  sensors  developed  under 
the  program  would  also  be  tested  for  long-range  capability.  Initial  analysis  of  the  long- 
range  experiments  indicates  that  the  target  should  be  elevated  to  avoid  the  extreme 
atmospheric  conditions  present  for  a  long  horizontal  path.  Two  options  for  target 
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location  exist:  the  target  can  be  located  on  a  distant  tower,  or  an  airborne  target  can  be 
used.  Certainly  the  airborne  target  can  provide  a  path  that  is  closer  to  vertical,  and  thus 
encompasses  lower  turbulence  levels.  However,  tracking  of  an  airborne  target  can  be 
more  challenging.  Furthermore,  the  class  of  airborne  targets  is  limited  by 
characteristics  of  the  heterodyne  receiver  array.  Specifically,  the  receiver  bandwidth 
and  frequency  offset  are  obtained  using  acousto-optic  methods  that  have  insufficient 
bandwidth  to  record  the  high  frequencies  inherent  to  high-speed  targets  such  as  rockets, 
unless  they  are  traveling  approximately  perpendicular  to  the  line  of  sight.  (This  would 
not  be  a  limitation  of  an  eventually  deployed  system.)  This  feature  drives  the 
experiments  to  use  stationary  or  slowly  moving  targets.  Of  the  several  candidates  for 
stationary  airborne  targets  (tethered  balloons,  helicopters,  remotely  piloted  vehicles 
(RPV),  etc.),  our  preliminary  analysis  shows  that  a  good  candidate  is  an  RPV  in  the 
form  of  a  radio-controlled  helicopter  with  the  target  tethered  below  it.  These  RPVs  are 
capable  of  flying  at  altitudes  limited  to  first  order  by  the  ability  to  see  the  RPV.  We 
estimate  that  such  an  RPV  can  be  piloted  by  using  a  telescope  to  an  altitude  of  3  km. 
The  target,  which  can  have  a  weight  of  roughly  3  lb,  is  tethered  at  a  distance  of  up  to 
100  ft  below  the  RPV.  The  tether  would  be  designed  to  damp  vibrations  from  the  RPV 
engine.  For  greater  altitudes,  a  manned  helicopter  can  be  used,  however,  the  expense  is 
considerably  greater  than  for  an  RPV.  Other  options  that  we  have  considered  include 
tethered  and  non-tethered  balloons;  however,  these  approaches  are  not  as  controllable. 
Sounding  rockets  may  be  used,  but  their  cost  limits  the  number  of  experiments  that  can 
be  performed. 

One  desirable  feature  associated  with  using  a  vertical  path  (airborne  target)  is  that 
the  level  of  turbulence  is  significantly  lower  than  for  a  horizontal  path.  Thus  the 
atmospherics,  the  anisoplanatic  patch  size  in  particular,  may  be  less  of  a  problem  for 
airborne  experiments.  It  is  anticipated,  however,  that  tracking  would  be  a  significant 
complication  for  the  airborne  case. 

Because  of  the  long  ranges  involved,  these  experiments  would  be  made  difficult  by 
the  lower  signal  levels.  These  signal  levels  can,  however,  be  offset  by  using  larger 
receiver  optics.  If  the  receiver  aperture  is  scaled  with  the  speckle  size,  the  signal  level 
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is  independent  of  range.  Thus  for  long  ranges  we  would  anticipate  replacement  of  the 
fore-optics  for  the  various  imaging  sensors. 

2.3.4.2  Experiment  Execution 

We  anticipate  that  these  long-range  experiments  would  make  significant  use  of 
hardware,  software  and  capabilities  developed  for  previous  experiments.  The  effort 
would  thus  concentrate  on  customizing  the  hardware,  such  as  mounting  the  sensors  on 
tracking  mounts  and  modifying  collection  optics,  and  planning  and  execution  of  the 
experiments. 
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3  TRADE-OFF  BETWEEN  ARRAY  SPARSITY  AND  NUMBER  OF 

REALIZATIONS 


The  quality  of  a  reconstructed  image  depends  on  numerous  factors,  including  signal- 
to-noise  ratio  (SNR),  aperture  diameter,  aperture  fill,  distance  to  target,  target 
structure,  wavelength,  phase  errors,  scattered  light,  etc.  In  this  section  we  show  digital 
simulations  demonstrating  the  trade-offs  for  (a)  different  amounts  of  aperture  fill,  (b) 
different  numbers  of  speckle  realizations,  and  (c)  different  types  of  objects.  We  do  this 
both  for  noise-free  diffraction-limited  coherent  imaging  and  for  noise-free  imaging 
correlography. 

For  a  review  of  coherent  imaging  and  imaging  correlography,  see  Appendix  A. 

3.1  DIFFRACTION-LIMITED  COHERENT  IMAGING 

Images  from  coherent  systems  also  suffer  from  a  mottled  appearance  that  we  call 
speckle,  which  decreases  image  intelligibility.  We  can  decrease  the  image  speckle  by 
collecting  several  frames  of  data,  each  with  the  object  having  translated  or  rotated 
slightly,  which  would  result  in  different  realizations  of  the  speckle  pattern.  By 
averaging  the  intensities  of  these  speckled  images  we  arrive  at  an  image  with  reduced 
speckle.  Averaging  over  a  large  number  of  speckle  realizations  results  in  an  image  that 
is  equivalent  to  a  speckle-free  image  of  the  target. 

For  the  case  of  imaging  with  sparse,  distributed  arrays  of  detectors,  the  details  of 
the  shape  of  the  impulse  response  play  a  major  role  in  image  quality.  The  impulse 
response  from  a  filled  array  typically  has  a  bright,  narrow  mainlobe  that  contains  most 
of  the  energy  and  small  sidelobes.  For  a  sparse  array  having  a  length  comparable  to  the 
length  of  a  given  filled  array,  the  impulse  response  would  be  similarly  narrow,  but  it 
would  contain  only  a  small  fraction  of  the  energy,  and  most  of  the  energy  would  be  in 
the  sidelobes.  The  result  is  typically  an  image  having  similar  resolution  as  that  from  a 
filled  array,  but  having  much  poorer  contrast  and  many  artifacts.  Equivalently,  the 
optical  transfer  function  (OTF),  or  its  magnitude,  the  modulation  transfer  function 
(MTF)  drops  off  quickly  at  the  lower  spatial  frequencies  as  compared  with  the  MTF 
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for  a  filled  array,  which  typically  has  approximately  a  triangular  shape.  For  reducing 
speckle  artifacts,  for  distributed  arrays  we  need  a  greater  number  of  realizations  than 
for  a  filled  array  from  the  simple  point  of  view  that  we  are  measuring  fewer  speckles 
at  any  given  time.  Similarly,  we  are  measuring  a  proportionally  smaller  number  of 
photons. 

The  theoretical  signal-to-noise  ratio  for  an  average  of  K  speckle  images  is,  in  the 
regime  of  low  sensor  and  photon  noise  (i.e.,  dominated  by  coherent  speckle  effects) 
[3.1] 

SNR(u)  =  ||i(u)|  VK  SBP  IOTF(u)|  (3.1) 

where  u  is  a  2-D  spatial-frequency  coordinate,  ||i(u)|  is  the  magnitude  of  the  Fourier 
transform  of  the  incoherent  object  normalized  to  be  unity  at  the  origin,  SBP  is  the 
space-bandwidth  product  of  the  data  and  OTF(u)  is  the  optical  transfer  function  (the 
autocorrelation  of  the  aperture  when  normalized  to  unity  at  the  origin).  SBP  equals  Ns, 

the  number  of  subapertures,  if  the  system  does  not  oversample  the  speckle  pattern. 

2 

Suppose  we  have  a  square  filled  array  with  L  subapertures  on  a  side,  or  a  total  of  L 
subapertures.  Then  SBP  =  L2,  and  at  the  middle  spatial  frequencies  OTF  =  (1/2),  so 
SBP  OTF  =  L2/2.  For  an  efficient  sparse  array,  however,  we  might  have  Ns  =  3(L  -  1) 
subapertures  (see  the  section  on  array  design).  For  a  nonredundant  sparse  array, 
OTF(u)  =  1/NS.  For  a  typical  array  pattern  we  would  have  OTF(u)  =  1/NS  at  many  of 
the  middle  spatial  frequencies.  Then  for  the  sparse  array  SBP  OTF  =  Ns/Ns  =  1.  Hence 
the  factor  V  SBP  OTF  in  the  expression  for  SNR  is  L/V2  times  greater  for  a  filled 
array  than  for  a  minimally-filled  sparse  array.  Hence  to  make  up  for  this  we  would 
need  the  number,  K,  of  speckle  realizations  to  be  L/^2  times  greater  for  the  sparse 
array.  We  could  also  argue  that,  at  the  spatial  frequency  components  outside  the  DC 
component,  the  OTF  of  an  efficient  sparse  array  has  an  average  value  of  the  OTF  of 
about  2/Ns.  Hence  we  may  need  only  L/(2yj2)  *  L/2.8  times  as  many  speckle 
realizations  for  the  sparse  array.  Alternatively,  we  can  build  sparse  array  patterns  that 
are  not  minimally  redundant  but  purposely  have  some  redundancy.  For  example,  with 
4L  -  2  subapertures  we  can  guarantee  that  each  spatial  frequency  has  an  OTF  of  at  least 
2.  Then  the  factor  SBP  OTF  =  Ns2/Ns  =  2,  and  we  would  need  only  half  as  many 
realizations  as  with  the  minimally  redundant  pattern  of  3L  -  3  subapertures. 
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After  averaging  image  intensities  to  form  a  speckle-reduced  image,  we  may  apply  a 
Fourier-domain  filter  that  boosts  the  higher  spatial  frequencies  so  that,  in  the  case  of 

I 

high  signal-to-noise  ratio  (SNR)  it  can  regain  the  character  of  an  incoherent  image 
through  a  filled  aperture.  This  factor  makes  sparse  arrays  particularly  appealing:  with 
enough  averaging  and  with  filtering  we  can  overcome  both  the  speckle  effects  and 
reduce  the  large  sidelobes  to  retrieve  a  speckle-free  image  that  is  as  good  as  one  that 
would  have  come  from  a  filled  aperture  costing  several  times  as  much. 

By  computer  simulations  we  demonstrated  the  trade-off  of  aperture  fill  versus 
number  of  speckle  realizations  for  different  classes  of  objects.  The  top  row  of  Figure 
3-1  shows  coherent  images  as  seen  through  a  filled  aperture  in  the  shape  of  a  circle. 
The  bottom  row  shows  images  of  the  same  object  as  seen  through  a  sparse  array 
consisting  of  21  small  circles,  arranged  in  an  annulus  with  outer  diameter  equal  to  the 
diameter  of  the  filled-aperture  circle.  The  object  consists  of  a  rectangle  of  size  2x3 
pixels  plus  four  additional  points,  each  1x1  pixel  in  size.  This  could  represent  a 
resolved  post-boost  vehicle  with  nearby  unresolved  re-entrant  vehicles  or  decoys.  We 
consider  this  “four-points”  object  to  be  a  fairly  simple  object.  An  ideal  (diffraction- 
limited  and  noise-free)  incoherent  image  of  the  object  is  shown  in  the  upper  right  of 
Figure  3-1.  The  slight  fuzziness  is  due  to  the  impulse  response  of  the  circular  aperture. 
A  single  coherent  image  of  the  object  is  shown  in  the  second  column  of  the  first  row. 
The  image  is  complex  valued,  and  the  figure  shows  the  intensity  of  the  image.  Note  that 
the  relative  brightnesses  of  the  points  is  distorted  and  the  2x3-pixel  block  at  the  top  is 
broken  up.  The  coherent  nature  of  the  imaging  process  produces  laser  speckles  in  the 
image  that  cause  these  effects.  Each  speckle  realization  would  break  up  and  distort  the 
image  intensity  in  a  different  way.  By  noncoherently  averaging  (i.e.,  averaging  the 
intensities  of)  images  with  several  different  speckle  realizations,  we  can  reduce  the 
speckle  artifacts  to  arrive  at  an  image  that  more  closely  resembles  the  incoherent 
image.  The  third  and  fourth  columns  show  the  effect  of  averaging  ten  and  one  hundred 
speckle  realizations,  respectively.  Even  with  as  few  as  ten  speckle  realizations  the  image 
looks  like  a  somewhat  noisy  version  of  the  ideal  incoherent  image  for  this  simple  object 
through  the  filled  aperture.  The  average  of  one  hundred  realizations  is  visually  nearly 
perfect. 
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Figure  3-1.  Imaging  Four- Points  Object  Through  a  Spai 


30 


Average  Ideal  Incoherent 
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Aperture  of  21  Circles. 


The  first  four  columns  of  the  second  row  of  Figure  3-1  show  the  same  thing  as  the 
first  row,  but  for  the  case  of  the  sparse  aperture  consisting  of  21  small  circles.  The 
single  coherent  image  is  almost  indecipherable.  The  average  of  ten  realizations  starts  to 
bring  out  some  of  the  features.  All  four  parts  to  the  object  finally  become  discernible 
after  averaging  100  realizations.  The  ring-like  artifact  that  is  especially  evident  around 
the  2x3  block  part  of  the  image  is  a  sidelobe  due  to  the  sparse  aperture.  This  artifact 
can  be  greatly  reduced  by  Wiener  filtering  [3.2],  the  result  of  which  is  shown  in  the 
rightmost  column  of  the  second  row.  In  the  process  of  filtering  the  sidelobes  are 
greatly  reduced,  but  at  the  cost  of  slightly  coarser  resolution  and  a  slight  enhancement 
of  the  background  noise.  The  quality  of  the  Wiener-filtered  average  of  100  images 
through  the  sparse  aperture  is  slightly  better  than  the  average  of  only  ten  realizations 
from  the  filled  aperture.  This  is  consistent  with  the  fact  that  there  is  about  five  times  as 
much  area  covered  by  the  filled  aperture  as  by  the  sparse  aperture. 

We  can  appreciate  the  effect  of  the  sparse  aperture  better  by  looking  at  Figure  3-2. 
The  first  column  shows  the  original  object  without  the  effects  of  speckle  or  low-pass 
filtering  by  the  aperture.  The  second  column  shows  the  speckle  pattern  in  the  aperture 
plane  caused  by  a  coherent  realization  of  the  object.  This  and  the  next  two  columns  are 
shown  demagnified  by  a  factor  of  two  compared  with  the  aperture  functions  shown  in 
Figure  3-1.  The  21 -subaperture  sparse  array  collects  the  speckles  shown  in  the  third 
column.  The  aperture  pattern  itself  constitutes  the  coherent  transfer  function  that 
defines  the  impulse  response  for  coherent  imaging.  Note  that  the  patterns  evident  in  the 
speckle  pattern  that  would  have  been  evident  in  a  filled-aperture  collection  are  less 
evident  in  the  sparse  array  collection.  The  right- most  column  shows  the  optical  transfer 
function  (OTF)  of  the  sparse  array  pattern.  This  is  the  transfer  function  for  the 
noncoherently  averaged  image.  The  central  part  was  overexposed  lest  the  surrounding 
parts  be  too  dark  to  see  clearly.  This  demonstrates  the  property  of  sparse  aperture  that 
the  high  spatial  frequencies  have  little  energy  compared  with  the  low  spatial 
frequencies.  Furthermore,  the  detailed  structure  of  the  OTF  causes  distortions  of  the 
objects  Fourier  transform  that  cause  a  distortion  in  the  image.  The  Wiener  filtering 
step  is  essentially  dividing  by  the  OTF  in  order  to  compensate  for  its  distortion.  The 
Wiener  filter  also  includes  a  noise  term  that  prevents  us  from  unduly  enhancing  the 
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noise  by  dividing  by  too  small  a  number  where  the  OTF  is  small.  This  is  described 
further  at  the  end  of  this  section. 

Figure  3-3  shows  similar  images  for  a  second  object,  consisting  of  a  cylinder  (post¬ 
boost  vehicle)  with  five  cones  (re-entrant  vehicles)  attached  plus  a  sixth  cone  detached 
from  the  cylinder.  This  object,  which  we  call  “bus+rv”  can  be  thought  of  as  a  much 
finer-resolution  image  of  the  central  portion  of  the  four-points  object,  for  which  the 
cylinder  is  now  well  resolved,  allowing  us  to  count  the  number  of  cones  attached  and  to 
discern  their  shape.  Note  that,  from  a  single  coherent  image  (top  row,  second  column) 
from  a  filled  aperture,  the  speckle  artifacts  make  it  impossible  to  identify  the  attached 
cones.  Only  the  overall  shape  can  be  easily  discerned.  The  noncoherent  average  of  one 
hundred  image  realizations,  shown  in  the  third  column,  yields  a  good-quality  image  that 
is  a  somewhat  noisy  version  of  the  ideal  incoherent  image.  Recall,  however,  that  for  the 
four-points  object  the  average  of  100  realizations  yielded  an  almost  perfect  image.  The 
bottom  row  shows  the  results  through  the  same  sparse  aperture  as  before.  The  single 
coherent  image  shown  in  the  second  column  is  almost  totally  useless — the  existence  of 
the  detached  cone  is  not  even  evident.  The  noncoherent  average  of  one  hundred  images, 
shown  in  the  third  column,  contains  some  information,  and  it  is  improved  considerably 
by  Wiener  filtering,  as  shown  in  the  fourth  column.  We  can  clearly  see  the  separated 
cone  (although  it  may  be  difficult  to  tell  that  it  is  a  cone),  and  the  cones  attached  to  the 
cylinder  can  be  seen,  but  counting  them  is  difficult.  Hence  we  see  that,  for  the  same 
aperture  and  the  same  number  of  image  realizations,  the  image  of  the  four-points  object 
has  considerably  higher  quality  than  the  image  of  the  bus+rv.  We  can  understand  this 
from  the  spatial  power  spectra  of  the  two  objects,  which  are  shown  in  Figure  3-4.  The 
spatial  power  spectrum  (and  its  square  root,  which  is  the  magnitude  of  the  Fourier 
transform  of  the  underlying  incoherent  object)  of  the  bus+rv,  which  is  an  extended 
object,  falls  off  much  more  rapidly  than  that  of  the  four-points;  hence  the  signal-to- 
noise  ratio  (SNR)  at  the  higher  spatial  frequencies  is  much  lower  for  the  bus+rv  object 
since  there  is  just  less  signal  there. 
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Figure  3-2.  Transfer  Functions  for  Imaging  Through  a  Sparse  Aperture  of  21  Circles. 
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Filled  Aperture  Coherent  Image  images  Incoherent  Image 


Sparse  Aperture  Coherent  Image 


Average 
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Figure  3-3.  Imaging  Bus+Rv  Object  Through  a  Sparse  Aperture  of  21  Circles. 


34 


Four  points 

Satellite 

Bus/RV 


The  examples  above  with  the  aperture  consisting  of  21  small  circles  are  for  the  case 
of  each  small  circle  having  several  detectors  across  it.  Hence  the  total  number  of 
detectors  is  several  times  21.  Suppose,  on  the  other  hand,  that  the  total  number  of 
detectors  is  a  similar  small  number,  only  22.  Figure  3-5  shows  the  result  for  a  small 
quasi-random  array  consisting  of  just  22  detectors.  Henshaw  [3.3]  designed  this  array 
configuration  using  a  genetic  algorithm  (see  Section  8).  Note  that  the  apparent 
resolution  is  many  times  coarser  than  for  the  21 -circles  aperture  because  of  the  great 
difference  in  numbers  of  detectors.  Figure  3-5  shows  the  bus+rv  and  4-points  images, 
at  a  lower  resolution,  as  seen  through  the  small-array  aperture.  Note  that  the 
noncoherent  average  of  100  realizations  of  the  bus+rv  image,  even  after  Wiener 
filtering,  shows  no  evidence  of  the  detached  cone.  Even  for  the  incoherent  image 
(which  is  equivalent  to  averaging  an  infinite  number  of  speckle  realizations),  shown  in 
the  right-most  column,  there  is  no  evidence  of  the  cones  on  the  cylinder.  These 
examples  demonstrate  a  fundamental  principal:  to  gather  an  image  with  a  fixed  array  of 
detectors,  we  need  a  large  number  of  detectors  to  arrive  at  an  image  with  a  large 
number  of  resolution  elements  across  it.  Hence  with,  say,  32  or  even  64  detectors  in  a 
fixed  pattern  we  can  hope  to  reconstruct  only  fairly  simple  images.  The  bottom  half  of 
Figure  3-5  shows  the  4-point  object  through  the  same  small  aperture;  in  this  case  we 
can  discern  the  shape  of  this  simple  object  with  a  reasonable  number  of  speckle 
realizations. 

If  we  allow  the  detectors  to  move  and/or  use  aperture  synthesis  by  observing  the 
speckles  as  they  flow  past  the  detectors,  then  extended  objects  can  be  successfully 
imaged  with  fewer  detectors,  but  with  a  much  higher  detector  read-out  rate. 
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Figure  3-5.  Imaging  Through  a  22-Element  Detector  Array; 

Top  row:  Extended  object  (“Bus/RV”),  simple  object  (“four  points”),  sparse  aperture, 
OTF.  Second  row:  noncoherent  average  of  coherent  images  of  Bus/RV  object.  Third 
row:  Wiener  filtered  images  of  Bus/RV.  Fourth  row:  average  of  images  of  four-points 
object.  Fifth  row:  Wiener  filtered  images  of  four-points  object. 
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3.2  IMAGING  CORRELOGRAPHY 


Some  of  the  advantages  of  imaging  correlography  are  that  it  allows  the 
reconstruction  of  a  fine-resolution  image  (a)  despite  atmospheric  turbulence, 
misaligned  array  elements,  and  other  sources  of  phase  errors,  and  (b)  with  relatively 
inexpensive  optical  hardware  and  detectors.  It  computes  the  average  autocovariance  of 
multiple  aperture-plane  speckle  intensity  realizations  to  estimate  the  magnitude  of  the 
Fourier  transform  of  the  incoherent  image  and  employs  a  phase  retrieval  algorithm  to 
reconstruct  the  image.  See  Appendix  A  for  more  details.  In  this  section  we  show  the 
effect  of  a  sparse  aperture  on  imaging  correlography. 

An  estimate  of  the  Fourier  magnitude  of  the  4-points  object,  determined  from 
correlography  using  100  speckle  realizations  through  the  21 -circle  aperture  (shown 
earlier),  is  shown  on  the  left  in  Figure  3-6.  The  higher  spatial  frequencies  are 
suppressed  because  of  the  fall-off  in  the  OTF  as  seen  in  Figure  3-2.  A  Wiener-filtered 
version  of  that  estimate,  in  which  the  suppressed  spatial  frequencies  are  enhanced,  is 
shown  in  the  second  column  of  Figure  3-6.  For  comparison,  the  third  column  shows  the 
true  Fourier  magnitude  of  the  object.  From  this  we  can  see  that  the  Wiener-filtered 
estimate  is  somewhat  depressed  at  the  lower  spatial  frequencies,  an  effect  we  have  not 
fully  explained.  The  right-most  column  of  the  Figure  shows  an  image  reconstructed 
from  the  Fourier  magnitude  estimate  (see  Appendix  A  for  an  explanation).  This 
reconstruction  is  a  faithful  representation  of  the  ideal  image,  demonstrating  that 
imaging  correlography  can  work  with  sparse  arrays  of  detectors. 

Figure  3-7  shows  the  effect  of  the  number  of  realizations  on  the  quality  of  the 
reconstructed  image.  The  top  row  shows  the  Fourier  magnitudes  (left-to-right):  of  the 
true  object,  estimated  from  one  hundred  speckle  realizations,  from  twenty  speckle 
realizations,  and  from  five  speckle  realizations.  The  bottom  row  shows  the 
corresponding  images  reconstructed  by  phase  retrieval.  For  twenty  or  more  speckle 
realizations  the  basic  structure  of  the  object  is  preserved  whereas  for  five  realizations 
the  image  is  badly  distorted. 
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Figure  3-6.  Fourier  Magnitude  and  Reconstructed  Image  from  Correlography. 
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Figure  3-7.  Image  from  Correlography  for  Different  Numbers  of  Realizations. 
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As  discussed  in  Appendix  A,  the  signal-to-noise  ratio  for  correlography  is  given  by 


SNR(u,  v)  = 


•\j N  ns  S(u,  v)  l|xl2 


where  nd  =  the  average  number  of  photons  per  detector 

M  =  the  number  of  detectors  (pixels)  per  aperture-plane  speckle  area 
M  nd  =  the  average  number  of  photons  per  speckle 
N  =  the  number  of  speckle  realizations 
ns  =  the  number  of  speckles  within  the  area  of  the  aperture 
f-t(u»  v)  =  Fj(u,  vj/F^O,  0),  and 

Fj(u,  v)  =  the  Fourier  transform  of  the  underlying  incoherent  object 


(3.2) 


9  9 

For  the  higher  spatial  frequencies  where  |fi|  is  small,  the  SNR  of  |p|  is  proportional  to 
N,  to  ||j.|  ,  and  to  the  square  root  of  the  area  of  the  aperture  (which  is  proportional  to 
ns).  However,  we  have  an  additional  proportionality  factor  of  S(u,  v),  the  OTF,  which 
also  depends  on  the  aperture  area  and  on  the  shape  of  the  aperture  pattern.  For 
example,  for  a  sparse  nonredundant  aperture  in  which  each  detector  sees  a  different 
speckle,  nsS(u,  v)  =  1  for  (u,  v)  *  (0,  0).  However,  for  a  sparse  redundant  aperture 
nsS(u,  v)  can  be  very  large.  In  general,  the  number  of  independent  Fourier  components 
sampled,  including  their  redundancy  but  ignoring  the  zero-frequency  component,  is 
Nd(Nd  -  l)/2,  where  Nd  is  the  number  of  detectors.  Hence  by  adding  just  one  more 
detector  to  a  given  number,  Nd,  we  add  (Nd  +l)(Nd)/2  -  Nd(Nd  -  l)/2  =  Nd  additional 
measurements.  The  quadratic  dependence  of  the  number  of  Fourier  measurements  on 
the  number  of  detectors  strongly  suggests  using  as  many  detectors  as  possible. 


3.3  SIMULATION  DETAILS 


Some  details  of  the  simulation  are  as  follows.  We  use  the  following  notation. 
fj(x)  =  incoherent  object 
Fj(u)  =  FT  of  incoherent  object 

9 

|Fj(u)|  =  power  spectrum  of  incoherent  object 
c(x)  =  coherent  impulse  response  (IPR) 
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s(x)  =  point-spread  function  (PSF)  =  |c(x)| 

A(u)  =  aperture  function  =  coherent  transfer  function  (a  binary  mask) 

OTF(u)  =  optical  transfer  function  =  autocorrelation  of  A(u) 

1.  Calculate  c(x),  s(x),  OTF(u). 

1.0  Make  sure  A(u)  is  of  width  <N/2  embedded  in  an  NxN  array 

1.1  c(x)  =  FT[A(u)] 

[if  you  want  to  use  c(x)  itself  later,  first  shift  A(u)  so  its  DC  is  in  u.l.  comer] 

1.2  s(x)  =  |c(x)|2  [imaginary  part  =  0] 

1.3  OTF^onnau^Cu)  =  Fr'[s(x)]  [Note:  FT1  includes  x  N'2  factor  ] 

1.4  OTF(u)  =  OTFunnormaijzetj(u)/OTFunnormajjze^(0) 

2.  Create  an  incoherent  image,  gr(x)  =  s(x)  *  fj(x),  from  f^x). 

2.1  fIc(x)  =  fj(x),  [imaginary  part  =  0] 

2.2  Fj(u)  =  FT[fIc(x)] 

2.3  Gr(u)  =  Fr(u)  OTF(u) 

2.4  s(x)*fI(x)  =  FT-1[GI(u)] 

2.5  Take  real  part  of  s(x)  *  fj(x)  [imaginary  part  should  «  real  part] 

3.  Create  a  coherent  image,  gk(x),  from  fj(x). 

3.1  The  real  and  imaginary  parts  of  gk(x)  are  chosen  from  a  zero-mean  Gaussian 
random  variable  with  variance  (1/2)  fj(x).  To  do  this, 

3.1.1  Form  a  uniform  array  of  complex  Gaussian  random  numbers 

(independent  Gaussian  r.n.  for  both  real  and  imaginary  parts)  each  having 
unity  variance  _ 

3.1.2  Multiply  that  array  by  ^  fj(x)/2  . 

(Note:  then  E[|gk(x)|2]  =  2fj(x)/2  =  fj(x).) 

3.2  FT[»]  [•  =  result  from  previous  step] 

3.3  •  x  A(u)  [having  DC  of  A(u)  in  u.l.  comer] 

3.4gk(x)=FT-1H 

Get  multiple  realizations  by  changing  the  value  of  the  seed  in  the  random  number 
generator  that  produced  the  Gaussian  r.n.’s.  To  be  able  to  reproduce  results,  always  use 
a  user-specified  seed. 
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4.  Wiener  filter  the  averaged  image. 

Filtering  an  image  g(x)  by  the  filter  function  H(u)  is  accomplished  by 
4. 1  G(u)  =  FT[g(x)] 

4-2  Gfilt(u)  =  G(u)  H(u) 

4-3  gfilt(x)  =  FT*1[G(u)] 

The  Wiener-Helstrom  filter  is  given  by 

OTFs*(u)  WF(u)  OTFf(u)  OTF *(u)  OTFf(u) 

Hw'h(U)  =  |OTFs(u)|2  WF(u)  +  WN(u)  =  |OTFs(u)|2  +  WN(u)/WF(u) 

OTFs*(u)  OTFf(u) 

=  |OTFs(u)|2  +  pn  NSR(u)  (33) 

where  OTFs(u)  is  the  OTF  of  the  sparse  aperture,  OTFf(u)  is  the  OTF  of  the  filled 
aperture,  WF(u)  =  |Fj(u)|2  is  the  power  (or  energy)  spectrum  of  the  object,  WN(u)  is 
the  power  spectrum  of  the  noise,  and  NSR(u)  =  WN(u)/WF(u)  is  the  noise-to-signal 
ratio  [l/SNR(u)],  and  pn  is  a  constant.  We  include  the  factor  OTFffii)  to  give  the  final 
image  an  aperture  weighting  closer  to  that  of  a  normal  filled  aperture.  Unfortunately, 
we  do  not  usually  know  WF(u)  and  sometimes  do  not  know  WN(u).  One  approach  is  to 
replace  NSR(u)  with  a  constant,  and  adjust  the  constant  to  give  the  most  pleasing  image 
(trading  off  image  smoothing  with  noise  suppression  for  image  sharpening  with  noise 
enhancement).  Using  a  constant  for  NSR(u)  usually  enhances  the  noise  at  the  highest 
spatial  frequencies  too  much.  Since  many  images  have  power  spectra  that  can  be 
approximated  by  a  power  law,  and  since  noise  power  spectra  are  often  a  constant,  a 
logical  choice  for  NSR(u)  would  be  a|u|X  where  a  and  y  are  constants,  y  controlling 
how  fast  the  WF(u)  falls  off  with  u.  Values  of  y  in  the  range  of  0.5  to  4  are 
appropriate,  with  the  smaller  values  appropriate  for  objects  consisting  of  a  very  small 
number  of  points  and  larger  values  appropriate  for  very  extended,  low-contrast  objects. 
This  expression  for  NSR(u)  must  be  modified  so  that  it  does  not  blow  up  near  u  =  0. 
We  accomplished  this  by  choosing 

f  1  ,  for  | u I  <  uL 

NSR(u)=  .  lY  ,  .  . 

[  a|u|T  ,  for  |u|  >  uF 
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where  a|uL|Y  =  1,  or  a  =  l/juj^.  Here  |u|  is  the  distance  in  pixels  from  the  DC  (i.e.,  if 
its  coordinate  relative  to  DC  is  (ux,  uy),  then  |u|  =  v/  ux  +  uy  .  Therefore  by 
specifying  uL  and  y,  NSR(u)  is  specified.  Examples  of  SNR(u)  and  NSR(u)  are  given  in 
Figures  3-8  to  3-11. 


SNR(u)  for  gam  =  1  and  ul  =  2,4,8,16 


Figure  3-8.  Model  of  Signal-to-Noise  Ratio  Versus  Spatial  Frequency.  For 
l/(frequency)  spectrum  and  model  parameter  ul  =  2,  4,  8,  and  16  pixels  (left  to  right). 
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NSR(u)  for  gam  =  2  and  ul  =  2,4,8,16 


Figure  3-11.  Model  of  Noise-to-Signal  Ratio  Versus  Spatial  Frequency.  For 
1 /(frequency)2  spectrum  and  model  parameter  uL  =  2,  4,  8,  and  16  pixels  (left  to 
right). 

As  long  as  the  OTF  does  not  include  aberrations,  then  OTF  =  OTF*  is  real-valued 
and  nonnegative. 
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4  PHASING  OF  SPARSE  ARRAYS  OF  HETERODYNE  RECEIVERS 


4.1  INTRODUCTION 

If  we  measure  an  optical  field,  scattered  from  a  laser-illuminated  object,  with  an 
array  of  heterodyne  receivers  [4.1],  then  we  can  reconstruct  an  image  of  the  object  by 
performing  a  suitable  propagation  transformation  (typically  a  Fresnel  or  Fourier 
transform)  in  the  computer  [4.2].  However,  if  phase  errors  are  present  in  the 
measurements,  then  a  blurred  image  will  result.  Phase  errors  can  be  induced  by  many 
things,  including  atmospheric  turbulence  and  path-length  errors  between  the  field  from 
the  object  and  the  local  oscillator  within  the  receiver.  Approaches  to  phasing  the  array 
without  a  beacon  either  require  the  array  of  heterodyne  receivers  to  be  on  a  filled 
regular  grid  [4.3,  4.4],  or  to  have  a  special  pattern  or  alternatively  require  a  low- 
resolution  image  be  available  [4.5]. 

In  this  section  we  describe  two  new  algorithms  for  phasing  heterodyne  arrays  that 
work  well  for  sparse  distributed  arrays  (and  should  work  for  filled  arrays  as  well), 
require  no  low-resolution  image,  and  are  relatively  insensitive  to  noise.  One  is  based  on 
the  iterative  transform  algorithm  using  a  support  constraint  [4.6,  4.7]  and  the  other  on 
maximizing  image  sharpness  [4.8,  4.9], 

4.2  IMAGING  MODEL 

With  the  heterodyne  receivers  we  measure  the  optical  fields  (ignoring  measurement 
noise) 

Gdk(u)  =  Fk(u)  exp[i<))e(u)]  ,  (4.1) 

for  k  =  1,  ...  ,  K,  where  Fk  is  the  ideal  complex  field  (without  phase  errors)  for 
speckle  realization  k,  4>e  is  the  phase  error,  and  u  is  a  2-D  coordinate  (the  pixel  index) 
in  the  measurement  plane.  Note  that  we  assume  the  phase  error  is  the  same  for  all  K 
realizations.  Different  speckle  realizations  are  obtained  if  the  object  rotates  slightly  or 
translates.  If  we  estimate  the  phase  error  as  <()(u),  then  our  estimates  of  the  optical  fields 
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Gk(u)  =  Gdk(u)  exp[-i<j)(u)]  (4.2) 

our  coherent  images  are 

gk(x)  =  J_1[Gk(u)]  =  N_2£Gk(u)  exp(i27tu-x/N)  ,  (4.3) 

u 

where  ux  is  a  vector  dot  (inner)  product,  and  our  speckle-reduced  averaged  image  is 

gl(x)  =  K-^lgkWI2  ,  (4-4) 

k 


where  summation  over  k  is  for  k  =  1,  ...,  K. 

Our  goal  is  to  estimate  the  image  [and,  equivalently  <j)e(u)  and  the  set  of  Fk(u)]  from 
the  set  of  Gdk(u). 

4.3  GENERALIZED  ITERATIVE  TRANSFORM  ALGORITHM 

Our  first  method  to  estimate  the  phase  error  is  a  generalization  of  the  iterative 
transform  algorithm  [4.6],  It  iterates  back  and  forth  between  the  image  and  Fourier 
domains,  applying  applicable  constraints  in  each  domain.  Only  here  we  are  performing 
the  iterations  simultaneously  on  all  K  image  realizations. 

Our  image-domain  constraint  is  a  support  constraint.  We  have  several  ways  to 
arrive  at  a  support  constraint. 

We  can  compute  the  averaged  autocorrelation  function,  and  use  our  methods  for 
estimating  image  supports  from  autocorrelation  supports  [4.10,  4.11].  Let 

rk(x)  =  ^-1[|Gdk(u)|2]  (4.5) 

be  the  autocorrelation  function  of  gdk(x)-  It  is  identical  to  the  unaberrated 
autocorrelation  since  it  depends  only  on  the  magnitude  of  Gdk(u)  and  not  on  its  phase. 
The  averaged  autocorrelation  function  is 

rI(x)  =  K-1Xlrk(x)|2.  (4.6) 

k 
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Except  for  an  impulse-response-like  term  at  the  origin,  rT(x)  will  approach  the  ideal 
incoherent  autocorrelation  for  large  values  of  K.  Note  that  this  computation  is 
essentially  the  same  one  we  use  for  imaging  laser  correlography  [4.12],  as  we  are  using 
only  the  intensity  of  the  measured  optical  fields.  It  is  important  to  use  the  averaged 
autocorrelation  rather  than  the  autocorrelation  from  an  individual  realization  since  any 
individual  autocorrelation  is  speckled  and  will  have  multiple  drop-outs  over  its  entire 
area,  making  its  support  difficult  to  estimate  [4.13]. 

We  may  also  estimate  the  object  support  by  using  a  low-resolution  image.  For 
example,  the  averaged  image  gi(x)  (the  current  blurred  estimate)  can  be  a  basis  for  a 
support  constraint.  Since  during  the  middle  of  the  reconstruction  sequence  gj(x)  will  be 
a  poor  approximation  to  the  true  image,  we  can  only  use  it  loosely  as  a  support 
constraint.  We  should  therefore  convolve  gj(x)  with  a  smoothing  function,  w(x),  giving 
s(x)  =  gj(x)  *  w(x),  which  also  has  the  advantage  of  being  a  relatively  smooth  function. 
We  can  use  this  as  a  gray-level  support  constraint  that  emphasizes  the  image  where  the 
probability  of  it  being  non-zero  is  highest  and  suppresses  it  where  the  probability  of  it 
being  nonzero  is  lower.  We  can  also  threshold  this  function  to  form  a  binary  mask 
representing  the  support. 

However  we  arrive  at  the  support  constraint,  we  will  denote  it  by  the  function  s(x), 
which  is  zero  outside  the  support  of  the  object.  To  apply  the  support  constraint,  we 
multiply  each  of  the  images  by  s(x): 

gsk(x)  =  s(x)  gk(x)  .  (4.7) 

After  transforming  these  support-constrained  images  back  to  the  Fourier  domain  we 
have  the  new  Fourier  estimates  Gsk(u)  =  7Tgsk(x)]- 

Our  Fourier  constraint  is  the  collection  of  magnitudes  of  the  measured  optical  fields, 
|Gdk(u)|.  The  phases  of  Gdk(u)  are  assumed  to  be  corrupted.  Furthermore,  in  the 
Fourier  domain  we  also  have  the  constraint  that  the  phase  error  is  common  to  all  K 
realizations.  That  is,  any  phase  changes  we  make  should  be  the  same  for  all  K 
realizations.  This  constraint  may  be  especially  strong  since  it  amounts  to  having  K  sets 
of  equations  to  solve,  but  only  one  set  of  phases  for  which  to  solve.  That  is,  compared 
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to  the  usual  phase  retrieval  problem,  this  one  may  be  much  more  heavily 
overdetermined,  which  may  make  it  more  robust.  On  the  other  hand,  if  the  array  data 
is  very  sparse,  then  we  may  be  at  a  disadvantage  compared  with  the  more  usual  filled- 
aperture  phase  retrieval  problems,  since  sparse  arrays  imply  large  sidelobes  that  would 
work  against  the  effectiveness  of  a  support  constraint.  The  step  of  subtracting  a 
common  phase  error  is,  in  effect,  using  a  Fourier  magnitude  constraint  since  the 
Fourier  magnitude  is  unchanged  for  each  k. 

The  new  Fourier  domain  estimates,  Gsk(u),  will  differ  from  the  previous  estimates 
Gk(u)  both  in  magnitude  and  in  phase.  In  the  usual  manner  for  an  iterative  transform 
algorithm,  we  substitute  for  the  magnitudes  of  Gsk(u)  the  measured  magnitudes 
|Gdk(u)|.  For  the  phases,  we  want  to  find  the  phase  functions  to  change  all  the  Gk(u)  by 
such  that  their  phases  most  closely  match  those  of  the  Gsk(u),  in  order  to  exploit  our 
knowledge  that  the  phase  error  is  the  same  for  all  realizations.  One  way  to  do  this  is  a 
follows.  Estimate  the  phase  error  by  taking  the  average  of  the  difference  in  phase 
between  Gdk(u)  and  Gsk(u).  At  a  given  point  (u),  for  some  of  the  realizations  Gsk(u) 
will  be  nearer  to  zero,  and  therefore  the  phases  for  those  realizations  will  be  less 
reliable  than  for  other  realizations.  Therefore  we  would  do  well  to  compute  a  weighted 
average  of  the  phase  difference.  Such  a  weighted  average  would  be 

4>(u)  =  arg{£  Gdk(u)  Gs^(u)}  ,  (4.8) 

k 

where  arg{X}  is  the  phase  of  X. 

In  summary,  one  iteration  of  the  generalized  algorithm  consists  of  the  following 
steps: 

(1)  Calculate  the  Fourier  estimates:  Gk(u)  =  Gdk(u)  exp[-i(|)(u)]  (for  k  =  1  ...  K). 

(2)  Inverse  FFT  them  to  compute  the  images  gk(x). 

(3)  Apply  the  support  constraint:  gsk(x)  =  s(x)  gk(x). 

(4)  FFT  to  compute  the  Gsk(u). 

(5)  Compute  a  new  phase-error  estimate:  <|)(u)  =  arg{£  Gdk(u)  Gs=|(u)}. 

k 

For  the  initial  iteration  the  phase  error  estimate  is  zero. 
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Steps  (2)  to  (4)  can  also  be  performed  by  a  direct  convolution  in  the  measurement 
plane: 

GSk(u)  =  Gk(u)  *  S(u)  =  £  S(Ul)  Gk(u  -  Ul)  ,  (4.9) 

ui 

where  S(u)  =  J  [s(x)],  which  would  be  faster  to  compute  for  the  case  of  sparse  arrays. 
For  filled  arrays  it  is  ordinarily  fastest  to  compute  a  convolution  using  two  FFT’s  and  a 
multiply.  However,  for  sparse  arrays,  for  which  Gk(u)  is  zero  except  at  P  locations, 
and  we  only  care  about  Gsk(u)  at  those  same  P  locations,  then  the  direct  convolution  to 
compute  Gsk(u)  will  take  only  KP2  complex  multiplications  [P  products  for  each  of  P 
values  of  Gsk(u)],  which  will  be  more  efficient  than  the  FFT’s.  The  FFT  method,  on  the 
other  hand,  is  very  straightforward  to  compute,  whereas  the  direct  convolution  is 
trickier  if  we  are  to  compute  it  efficiendy. 

If  we  were  to  choose  for  the  support  constraint  s(x)  =  gT(x)  *  w(x),  then  S(u)  = 
Gj(u)  W(u),  and 

Gsk(u)  =  Gk(u)  *  [Gjfu)  W(u)]  =  £  G^UjlWfu!)  Gk(u  -  Ul)  .  (4.10) 

ui 

Since  w(x)  would  ordinarily  be  a  real-valued  symmetric  kernel  of  width  a  small 
number  of  pixels  (depending  on  how  close  we  think  we  are  to  a  solution),  W(u)  would 
be  a  smooth  real-valued  symmetric  function  that  would  play  the  role  of  a  weighting 
function  in  the  Fourier  domain. 

4.4  MAXIMIZING  SHARPNESS 

The  second  approach  is  to  find  the  phase-error  estimate  0(u)  which  maximize  the 
sharpness. 

Si  =  X[gi(x)]2  =  I[Xlgk(x)l2]2  ,  (4.11) 

X  X  k 


of  the  image  estimate. 
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We  will  accomplish  the  maximization  using  a  standard  gradient  search  technique  that 
employs  an  analytic  expression  for  the  gradient,  which  we  derive  as  follows: 


asi 

a<j>(u) =  2 


£gl(x)X 

x  k 


dlgkPQI2 

a<!>(u) 


(4.12) 


Using 

d|ty(u)  =  _iN_2  §k*(x)  Gk(u)  exp(i2tcu-x/N)  +  c.c.  ,  (4.13) 

where  c.c.  denotes  the  complex  conjugate  of  the  term  that  precedes  it,  and  rearranging 
the  order  of  summations  yields 

dS 

=  4N~2  Im{  XGk(u)  (£gk(x)  gi(x)  exp(-i2jiu-x/N))*  } 

=  4N~2  Im{  £Gk(u)  (^[gkW  glWl)*}  •  (4-14) 

k 

Sj  would  be  maximized  by  a  gradient  search  algorithm  employing  the  expressions 
above  for  S  j  and  its  gradient. 

Note  that  the  product  gk(x)  gj(x)  resembles  constraining  the  support  of  gk(x)  to  fit 
within  the  support  of  gj(x).  Hence  this  algorithm  is  somewhat  related  to  the  iterative 
transform  algorithm  and  other  phase-retrieval  gradient-search  algorithms  we  have 
employed.  We  give  a  more  detailed  analysis  of  this  point  later. 

If  we  are  dealing  with  sparse  apertures,  gj(x)  will  have  large  sidelobes,  and  we  may 
wish  to  filter  it  first  to  undo  the  severe  weighting  due  to  the  optical  transfer  function. 
Let  a  filtered  version  of  gj(x)  be 

giw(x)  =  ^“MGlw(u)l  =^_1{W(u)  Jtg^x)]} 

=  N-2  Xw(u')  exp(i27tu'-x/N)  Xgl(x')  exp(-i27tu'-x'/N)  ,(4.15) 
u'  x' 

where  W(u)  is  a  filter  function,  a  good  choice  for  which  would  be  a  Wiener-Hel strom 
filter.  W(u)  will  be  centro-symmetric,  making  gIW(u)  real-valued,  but  not  necessarily 
nonnegative.  Then  our  new  sharpness  metric  is 
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(4.16) 


Siw  -  X[giw«]2  • 

x 

Similar  to  the  approach  above,  we  take  the  derivative  of  with  respect  to  <t>(u), 
reverse  the  order  of  the  (four)  summations,  and,  using  the  fact  that  W(u)  is  Hermitian, 
we  obtain 

dS 

=  4N-2  Im{2Gk(u)  (5{gk(x)  r’W))*) 

=  4N“2  Im{XGk(u)  (Gk(u)  *  [G,W2])*)  ,  (4.17) 

k 

where  GIW(u)  =  W(u)  7  [gi(x)],  making  GIW(u)W(u)  =  W2(u)  7  [gi(x)].  This 
expression  is  the  same  as  that  for  Sj  given  above,  except  that  gj(x)  is  replaced  by  7 
_1[GiwW]  =  giw(x)  *  f-'m  =  gi(x)  *  y-l[W2]  ,  where  *  denotes  convolution.  S1W 
reduces  to  Sj  if  W(u)  =  1  everywhere. 

The  major  burden  of  computing  Sj  is  the  K  FFT’s  needed  to  compute  the  gk(x)’s. 
The  major  computational  burden  of  computing  the  gradient  of  Sj  with  respect  to  <}>(u) 
is  the  2K  FFT’s  needed  to  compute  gk(x)  and  yjgj(x)  gk(x)]. 

The  primary  burden  of  computing  is  (2K  +  2)  FFT’s,  which  is  only  slightly 
greater  than  that  of  computing  Sj. 

4.5  COMPARISON  WITH  MAXIMIZING  SHARPNESS 

In  this  section  we  show  that  the  error-reduction  [4.6]  version  of  our  generalized 
iterative  transform  algorithm  is  similar  in  effect  to  a  steepest-descent  gradient  search  to 
maximize  sharpness  with  a  fixed  step  size.  (However,  we  in  practice  use  a  more 
powerful  gradient  search  algorithm,  such  as  conjugate-gradient,  and  we  can  use  a  more 
powerful  version  of  the  iterative  transform  algorithm  as  well,  such  as  hybrid  input- 
output  [4.6]). 

In  the  section  above  we  derived  a  gradient  of  a  sharpness  merit  function,  given  as 
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H®=4N-2  Im{|;Gk(u)(7-(gk(x)  iT-l[G,w2]  })*}  ,  (4.18) 

for  the  weighted-Gj  case,  where  Gj(u)  is  replaced  by  GIW(u)  =  Gj(u)W(u).  The 
unweighted  case  would  be  the  same  with  W(u)  =  1  everywhere.  This  can  be  rewritten 

3S 

3^=  4N-2  Im{|,Gk(u)  (Gk(x)  *  [G,W2])*}  .  (4.19) 


Let 

Gtk(u)=  Gk(u)  *  [GjW2] .  (4.20) 

Note  that  G^u)  is  very  similar  to  Gsk(u)  —  it  is  identical  if  W(u)  is  a  binary  function. 
We  have 

dS 

3^™=4N-2  Im{XGk(u)Gft(u)|  .  (4.21) 

For  small  values  5,  Im{|A|  exp[i8(u)]}  *  |A|  5(u),  and  by  definition  arg{|A|  exp[i8(u)]} 
=  8(u),  or  Im{|A|  exp[i8(u)]}  =  |A|  arg{|A|  exp[iS(u)]}.  Therefore  in  the  equation 
above,  the  gradient  with  respect  to  the  phase  error  estimate  is  proportional  to  the 
average  phase  difference  between  Gk(u)  and  G^u).  Since  we  want  to  maximize  the 
sharpness  (climb  up  the  hill)  we  want  the  phase-error  estimate,  4>(u),  to  go  in  the 
direction  of  this  gradient.  So  if  the  phase  moves  in  that  direction  by  a  step  size,  b,  then 
the  new  phase  error  estimate  will  be 

<f>new(u)  =  <t>old(u)  +  b  arg{£Gk(u)  G£(u)} 

k 

*  $old(u)  +  b  arg  { £Gdk(u)  exp[-i<j>old(u)]  G£(u) }  .  (4.22) 

k 

where  we  have  included  the  magnitude  of  the  sum  in  the  step  size  b.  If  we  happen  to 
choose  a  step  size  b=l,  then  we  have,  for  maximizing  sharpness, 

<t)new(u)  “  arg  {  XGdk(u)  G^(u)}  •  (4-23) 

k 
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Comparing  this  with  our  new  phase  error  estimate  for  the  iterative  transform 
algorithm  (Step  5),  we  see  that  it  is  nearly  the  same,  since  G^(u)  ~  Gsk(u).  That  is,  our 
iterative  transform  algorithm  using  a  low-resolution  image  as  a  support  constraint  is 
approximately  the  same  as  a  sharpness-maximization  steepest-descent  gradient  search 
method  with  a  particular  fixed  step  size. 

Which  of  the  two  algorithms  is  better  is  difficult  to  predict.  We  know  from  previous 
work  that  the  error-reduction  (ER)  version  of  the  iterative  transform  algorithm  (ITA) 
is  equivalent  to  a  fixed- step  steepest  descent  algorithm,  and  it  is  inferior  to  both  the 
hybrid  input-output  (HIO)  version  of  the  ITA  and  to  a  conjugate-gradient  algorithm. 
The  algorithm  we  described  above  is  an  ER-like  algorithm.  Therefore  the  sharpness- 
maximizing  gradient-search  algorithm  (which  would  employ  a  sophisticated  gradient 
search)  would  be  superior  to  the  iterative  algorithm  given  above  both  because  the  type 
of  gradient  search  is  superior  and  because  it  is  not  restricted  to  a  specific  step  size.  On 
the  other  hand,  a  HIO  version  of  this  algorithm  might  work  well.  An  HIO  version  of 
the  iterative  algorithm  would  be  to  replace  Step  3,  gsk(x)  =  s(x)  gk(x),  with 

gsk(x)  =  s(x)  gk(x)  +  [1  -  s(x)]  [gskoid(x)  -  p  gk(x)] 

=  s(x)  [(1  +  p)  gk(x)  -  gskold(x)]  +  gskoid(x)  -  P  gk(x)  (4.24) 

where  gskoid(x)  is  gsk(x)  from  the  previous  iteration  and  we  assume  that  s(x)  is  binary 
(1  or  0),  although  this  might  also  work  for  s(x)  defined  continuously  between  0  and  1. 
A  Fourier-domain  convolution  version  of  the  HIO  algorithm,  which  would  speed  the 
algorithm  for  sparse  arrays,  would  be  as  follows  (replacing  Steps  2  to  4): 

Gsk(u)  =  S(u)  *  [(1  +  P)  Gk(u)  -  Gskold(u)]  +  Gskold(u)  -  p  Gk(u)  .  (4.25) 

4.6  IMAGE  RECONSTRUCTION  EXAMPLE 

Figure  4-1  shows  a  computer  simulation  example  of  phase-error  correction  using 
the  generalized  iterative-transform  algorithm  for  a  distributed  array  of  heterodyne 
receivers.  At  the  top  is  the  4-points  object,  consisting  of  a  block  of  2x3  pixels  and  three 
separated  points,  embedded  in  a  64x64  array.  The  sparse  aperture  we  used  for  this 
experiment  has  54  subapertures  within  a  16x16  area,  which  is  shown  in  the  upper  left, 
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magnified  by  a  factor  of  two  (only  the  central  32x32  pixels  in  its  64x64  array  are 
shown).  This  sub-aperture  pattern  was  designed  using  a  genetic-based  algorithm  by  P. 
Henshaw  [4.14].  The  average  of  ten  ideal  diffraction-limited  speckled  images  through 
this  aperture  has  significant  sidelobe  artifacts  making  it  difficult  to  discern  the  four 
points,  as  seen  in  the  next  part  of  the  figure.  The  Wiener- Helstrom  [4.15]  filtered 
version  of  the  averaged  image  (on  the  right)  has  few  of  these  artifacts,  and  from  it  one 
can  easily  discern  the  major  parts  of  the  image.  The  bottom  row  in  the  figure  shows  the 
results  when  we  added,  to  the  Fourier  transform,  random  phase  errors  that  were 
uniformly  distributed  Gaussian  random  numbers  with  standard  deviation  2ti  radians. 
The  raw  aberrated  image  (on  the  left)  has  no  information  content,  whereas  the  image 
corrected  by  the  iterative  transform  is  almost  identical  to  the  ideal  diffraction-limited 
image.  The  reconstructed  image  was  shifted  in  position  from  the  ideal,  as  is  usually  the 
case  for  phase-correction  algorithms.  (All  images  in  this  example  are  the  average  of  ten 
speckled  images.)  We  obtained  comparable  results  using  the  image  sharpness  algorithm. 
The  algorithms  usually  produced  images  with  better  quality  as  the  number  of 
realizations  increases.  It  was  unreliable  when  using  only  one  to  four  realizations.  Tests 
showed  that  the  algorithm  is  not  very  sensitive  to  noise. 

4.7  CONCLUSIONS 

We  have  derived  two  new  algorithms  for  phasing  an  array  of  heterodyne  receivers, 
and  have  demonstrated  that  they  can  work  well  for  sparse  distributed  arrays,  require  a 
modest  number  of  speckle  realizations,  and  are  robust  in  the  presence  of  measurement 
noise. 
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Sparse  array  Ideal  ifoage  Wiener-Filtered 

(no  phase  error)  idea1  Image 


Aberrated  Image  Corrected  Image  Wiener-Filtered 

(2/t  phase  error)  Corrected  Image 


Figure  4-1.  Computer  Simulation  of  Sparse  Heterodyne  Array  Phasing  by  the  Iterative 
Transform  Algorithm. 
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5  TURBULENCE  ISSUES  FOR  LASER-IMAGING  EXPERIMENTS 

5.1  INTRODUCTION 

Atmospheric  turbulence  plays  an  important  role  in  the  design  of  laser-imaging  ex¬ 
periments.  Such  experiments  could  include  the  investigation  of  heterodyne  methods, 
quasi-heterodyne  methods,  and  direct-detection  methods.  Phase  aberrations  caused 
by  atmospheric  turbulence  add  to  receiver  misalignments  in  heterodyne  detection  to 
cause  a  total  phase  error  that  can  cause  severe  blurring  of  the  imagery.  In  addition, 
anisoplanatic  effects  (phase  errors  being  different  for  different  points  in  the  object) 
are  not  accommodated  by  any  of  these  methods,  in  their  current  formulation.  It  is 
natural  to  consider  horizontal-path  geometries  in  order  to  simplify  the  experimen¬ 
tal  design.  However,  induced  phase  aberrations,  anisoplanatism,  and  scintillation 
are  much  more  severe  with  horizontal  geometries  than  with  vertical  geometries.  A 
candidate  experiment  design  might  set  the  range  at  R  =  1  km,  the  sparse  aperture 
diameter  at  D  =  25  cm  and  the  wavelength  at  A  =  1.0e-6m  to  get  a  target  resolution 
of 


Ax  =  RX/D 

(5.1) 

=  3.9  mm. 

(5.2) 

Phase-aberration  strength  (quantified  with  a  correlation  diameter)  and  anisopla¬ 
natism  (quantified  with  an  isoplanatic  patch  diameter)  should  be  evaluated  with 
respect  to  the  aperture  diameter  and  extent  of  the  target  (Ax  times  number  of  res¬ 
olution  elements  across  target),  respectively.  In  this  Section  we  investigate  the  ex¬ 
pected  phase-aberration  strength  and  anisoplanatism  for  a  horizontal-path  geometry 
at  ranges  of  1  and  10km. 

5.2  TURBULENCE  STRENGTH 

The  index  structure  parameter,  C%,  is  a  measure  of  turbulence  strength  that 
depends  upon  altitude,  time  of  day,  meteorological  conditions,  and  season.  Induced 


61 


®>ERIM 


Table  5-1.  Correlation  Diameter,  rD,  as  a  Function  of  C\  and  Range,  R. 


Cl  \  R 

rQ 

1km  10  km 

HU15 

IO-14 

10-13 

10-12 

33.3  cm  8.36  cm 
8.36  cm  2.10  cm 
2.10  cm  0.53  cm 
0.53  cm  0.13  cm 

phase  aberrations,  anisoplanatism,  and  scintillation  derive  from  path  integrals  involv¬ 
ing  (72,  under  the  Kolmogorov  turbulence  model.  Horizontal  imaging  implies  that  C'l 
is  independent  of  position  on  path,  simplifying  computation  of  these  path  integrals. 
Note  however  that  Cl  is  highly  variable,  so  that  several  computations  must  be  made. 
Figure  5-1  shows  the  variability  in  Cl  at  an  altitude  of  9  m  over  the  diurnal  cycle 
under  cloudless  desert  conditions.  This  figure  suggests  that  Cl  could  be  expected 
to  vary  over  4  orders  of  magnitude  at  a  given  site.  Fried  has  defined  the  correlation 
diameter  (or  seeing  parameter)  to  be  [5.1] 

(5-3) 

Using  this  expression  we  computed  the  correlation  diameter  as  a  function  of  Cl  and 
range.  These  values  are  given  in  Table  5-1.  We  see  from  Table  5-1  that  only  the  1  km 
range  and  the  most  favorable  turbulence  conditions  yield  a  correlation  diameter  that 
is  comparable  to  the  hypothetical  aperture  diameter,  D  —  25  cm. 

5.3  ANISOPLANATISM 

Another  important  factor  to  consider  in  an  experiment  design  is  the  amount  of 
anisoplanatism  induced  by  the  long-path  turbulence.  The  origin  of  the  anisopla¬ 
natism  is  the  fact  that  radiation  emanating  from  differing  points  in  the  field  of  view 
encounter  different  volumes  of  turbulence  in  route  to  the  aperture.  Fried  has  derived 
an  expression  for  anisoplanatism-induced  mean-square  phase  error  in  an  adaptively 
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Figure  5-1.  Cl  Measurements.  Measurements  are  taken  at  an  altitude  of  9m  under 
cloudless  desert  conditions  over  the  diurnal  cycle.  Statistical  variability  in  the  form 
of  ±1  standard  deviation  is  also  given  (from  Infrared  and  Electro-Optical  Systems 
Handbook,  SPIE  Optical  Engineering  Press,  Bellingham  WA,  Vol.  2,  p.  202). 
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corrected  system  [5.2]: 

5( r.a)  =  j*Cl(z)  ({r[l  -  (z/i?)]}5'3  +  M5/3 

-j|{r[l  -  (zIR)]?  +  2razc[l  -  ( z/R )]  +  (azf|s/6 
-||{r[l  -  (zlR)]f  -  2razc[l  -  (z/R)}  +  (az)2|5/6)  dz  (5.4) 

where  a  is  a  field-angle  vector  relative  to  the  beacon  and 

c  =  r  ■  a/(ra)  .  (5-5) 

The  function  S(r,a.)  is  interpreted  as  the  mean-square  phase  error  between  points  in 
the  pupil  separated  by  r  for  a  point  in  the  field-of-view  at  field-angle  a.  Equation  (5.4) 
applies  to  a  system  with  ideal  phase  correction  for  waves  generated  at  the  beacon. 
The  RMS  phase  difference  at  a  separation  of  25  cm  in  the  pupil  as  a  function  of  field 
position  and  at  a  range  R  =  1  km  with  c  =  1  is  given  in  Figure  5-2. 

The  isoplanatic  angle  has  been  defined  by  Fried  to  be 


Notice  that  this  expression  resembles  the  definition  of  r0  (Eq.  (5.3)),  although  there  is 
a  different  scale  factor  and  the  path  integral  emphasizes  turbulence  near  the  target. 
In  order  to  interpret  the  meaning  of  the  anisoplanatic  angle  it  is  useful  to  consider 
the  anisoplanatism  induced  mean-square  phase  error  in  two  important  limiting  cases. 


lim  5(r,a) 

r/a— >  oo 

lim  S(r,a ) 

r/a— >0 


(5.7) 

(5.8) 


We  define  the  isoplanatic  patch  diameter  as 


p0  =  Ra0  . 


(5.9) 


64 


®>ERIM 


Table  5-2.  Isoplanatic  Patch  Diameter  as  a  Function  of  Turbulence  Strength,  C„,  and 


Range,  R. 


Table  5-2  gives  the  isoplanatic  patch  diameter  as  a  function  of  turbulence  strength, 
C%,  and  range,  R.  Table  5-2  suggests  that  for  an  experiment  design  with  R  —  1  km  and 
very  good  turbulence  conditions,  we  might  expect  to  have  between  7  (for  C*  =  10~14) 
and  27  (for  C%  =  10“15)  resolution  elements  across  an  isoplanatic  patch  when  the 
resolution  is  Ax  =  3.9  mm. 

Consider  the  case  in  which  we  are  not  troubled  by  phase  aberrations  in  the  pupil 
but  are  only  constrained  to  have  the  field  of  view  confined  to  an  isoplanatic  patch. 
This  would  be  the  case  when  we  have  an  adaptive  optics  system  or  when  we  have  full 
confidence  in  a  phase-correction  algorithm.  In  this  case  we  have 

NAx  =  p0  (5.10) 

=  Ra0  .  (5.11) 


Using  Eq.  (5.1)  to  substitute  for  Ax  and  solving  for  D  we  get 


D 


NR\ 

(Ra0) 


(5.12) 


Let  N  =  31,  R  =  1km,  and  A  =  1.0 ^m.  Table  5-3  gives  the  required  aperture 
diameter  values  as  a  function  of  turbulence  strength,  C%. 
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Table  5-3.  Aperture  Diameter  as  a  Function  of  Turbulence  Strength,  C\.  Fixed 
parameters  include  N  =  31,  R  =  1  km,  and  A  -  1.0  /im. 


Cn 

Ra0 

D 

10-15 

105.0  mm 

10-14 

26.3  mm 

9 

10-13 

6.6  mm 

1' 

10~12 

1.7  mm 

18.82  m 

5.4  RECOMMENDATIONS  AND  OUTSTANDING  ISSUES 

These  calculations  have  been  made  using  plausible  values  for  C\.  We  recommend 
that  on-site  measurements  be  made  of  C\  to  help  guide  the  experiment  design.  We 
also  recommend  that  the  simulation  experiments  be  performed  to  investigate  the  ef¬ 
fects  of  anisoplanatism  on  phase  correction  algorithms.  If  necessary,  these  algorithms 
could  be  developed  to  explicitly  accommodate  anisoplanatic  effects.  It  appears  that 
a  horizontal  range  of  no  more  than  1  km  should  be  used  for  initial  experiments.  Data 
may  need  to  be  collected  at  the  diurnal  minima  for  C2  and  a  temperature  structure 
constant,  U2,  may  need  to  be  monitored  to  optimize  data-collection  times. 

There  are  several  outstanding  issues  needing  further  study.  The  effects  of  boundary- 
layer  turbulence  (due  to  proximity  to  the  ground)  have  not  been  considered.  On-site 
measurements  should  be  designed  to  statistically  characterize  C 2  as  a  function  of  time 
of  day  and  height  off  of  the  ground.  In  addition,  effects  of  scintillation  also  need  to 
be  characterized  for  any  experimental  design.  Finally,  an  extended-path  experiment 
might  be  designed  with  a  vertical  geometry.  Although  such  a  geometry  significantly 
reduces  problematic  atmospheric  effects,  target  design  becomes  more  challenging. 
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6  NEW  3-D  IMAGING  MODALITY  USING 
AN  OPACITY  CONSTRAINT 

6.1  INTRODUCTION 

In  this  section  we  present  a  candidate  three-dimensional  imaging  concept  that 
we  refer  to  as  Phase  Retrieval  with  an  Opacity  Constraint  for  LAser  IMaging  (PRO¬ 
CLAIM).  Three-dimensional  imaging  significantly  enhances  target  identification  over 
two-dimensional  imaging  modes. 

6.2  IMAGING  CONCEPT 

Marron  has  demonstrated  that  three-dimensional  imaging  can  be  accomplished 
by  sequentially  illuminating  an  object  with  different  laser  frequencies  and  measuring 
the  far-field  speckle  pattern  for  each  of  the  illuminating  frequencies  [6.1].  In-phase 
and  quadrature  measurements  provide  field  quantities  for  the  laser-speckle  patterns. 
Properly  formatted,  these  data  form  a  three-dimensional  Fourier-volume  (or  aperture) 
representation  of  the  illuminated  object.  A  simple  three-dimensional  DFT  can  be 
performed  to  get  a  three-dimensional  image  of  the  object.  This  three-dimensional 
lensless  imaging  concept  is  referred  to  as  Holographic  Laser  Radar  (HLR). 

The  collection  of  HLR  data  could  be  considerably  simplified  if,  instead  of  collecting 
field  measurements  for  the  speckle  images,  we  collect  intensity  measurements.  Inten¬ 
sity  measurements  are  straightforward  and  eliminate  the  need  for  interference  with 
a  reference  beam  with  precision  alignment,  tracking,  and  phase  stability.  A  speckle- 
intensity  data  set  would  provide  three-dimensional  Fourier-magnitude  information 
or,  equivalently,  the  three-dimensional  object  autocorrelation.  In  order  to  recover  a 
literal  three-dimensional  image,  a  phase-retrieval  algorithm  is  required.  Therefore, 
this  proposed  imaging  concept  trades  complexity  and  cost  in  hardware  for  increased 
computing. 

Phase  retrieval  requires  some  type  of  a  priori  information  about  the  object.  Two- 
dimensional  complex-valued  objects  have  been  recovered  using  phase  retrieval  with  a 
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support  constraint  [6.2],  although  this  is  a  challenging  problem.  It  is  well  known  that 
the  uniqueness  properties  of  two-dimensional  phase  retrieval  are  much  better  than 
for  the  one-dimensional  problem  [6.3].  We  conjecture  that  three-dimensional  phase 
retrieval  with  a  support  constraint  is  better  conditioned  than  its  two-dimensional 
counterpart. 

An  alternative  constraint  that  has  great  promise  in  the  three-dimensional  imaging 
case  is  an  opacity  constraint.  An  opaque  object  is  an  object  that  exhibits  only  surface 
scattering  and  no  volume  scattering  (over  volumes  that  extend  beyond  the  desired 
range  resolution).  The  reflectivity  function  for  an  opaque  object  is  confined  to  a  two- 
dimensional  manifold  embedded  in  a  three-dimensional  space.  An  opacity  constraint 
has  been  used  to  perform  superresolution  [6.4]. 

The  opacity  constraint  is  a  special  type  of  support  constraint.  It  is  a  “quality 
of  support”  constraint  -  the  actual  location  of  the  support  is  not  given,  although 
the  object  is  known  to  be  confined  to  a  two-dimensional  manifold.  This  constraint 
promises  to  be  very  powerful  since  it  greatly  reduces  the  class  of  feasible  objects 
from  which  to  choose  an  estimate.  Moreover,  there  are  many  imaging  applications  in 
which  the  objects  will  be  known  with  confidence  to  be  opaque.  For  example,  ballistic 
missiles  certainly  qualify  as  opaque  objects.  Most  objects  in  our  everyday  experience 
satisfy  the  opaque  condition.  The  constraint  is  invalid  for  objects  with  distributed 
volume  scatterers  such  as  translucent  or  fog-like  objects. 

The  use  of  an  opacity  constraint  in  conjunction  with  frequency-diverse  pupil- 
plane  speckle  measurements  to  reconstruct  a  three-dimensional  object  constitutes  the 
unconventional  imaging  concept  that  we  refer  to  as  Phase  Retrieval  with  an  Opacity 
Constraint  for  LAser  IMaging  (PROCLAIM).  In  the  following  sections  we  consider 
various  embodiments  of  this  imaging  concept. 


6.3  SYNTHETIC-APERTURE  DESIGN 

Consider  a  missile  or  a  space  object  that  is  moving  with  some  trajectory  overhead. 
When  it  is  illuminated  with  a  laser,  the  reflected  light  will  form  a  speckle  pattern  on 
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the  ground.  The  diameter  of  a  speckle  on  the  ground  is  given  by 


P  = 


RX 

D0bj  ’ 


(6.1) 


where  R  is  the  target  range,  A  is  the  wavelength,  and  D0bj  is  the  object  diameter 
in  the  cross-range  direction.  We  wish  to  detect  this  speckle  pattern.  The  spatial 
sampling  interval  required  to  satisfy  the  Nyquist  criterion  is 


Ax  <  p[  2 

(6.2) 

RX 

(6.3) 

~  2  D0bj 

For  an  object  with  D0bj  —  5  m  at  R  =  100  km  and  A  =  1  /jm  we  have  Ac  <  1  cm. 


One  way  to  collect  an  angle-angle  speckle  pattern  at  a  particular  frequency  is 
to  use  a  filled  array,  as  depicted  in  Figure  6-l(a).  Filled  arrays  can  be  expensive. 
An  alternative  is  to  recognize  that  the  object  will  be  rotating  and  the  speckle  will 
be  translating  on  the  ground,  which  suggests  the  use  of  a  synthetic  aperture.  The 
speckle  will  translate  on  the  ground  at  a  speed  given  by 


vs  —  2Rlo  (6-4) 

where  u>  is  the  angular  velocity  of  the  target.  If  a  linear  array  of  detectors  is  oriented 
perpendicular  to  the  direction  of  speckle  motion,  an  angle- angle  speckle  pattern  can 
be  built  up  synthetically  in  time,  as  depicted  in  Figure  6-1  (b).  The  hardware  require¬ 
ments  for  such  a  system  are  dramatically  reduced  over  a  filled  array  of  comparable 
diameter. 


6.4  SERIAL-FREQUENCY  ILLUMINATION 

The  goal  is  to  acquire  a  separate  array  of  speckle  measurements  for  each  of  several 
illuminating  laser  frequencies.  This  can  be  accomplished  by  modulating  the  illu¬ 
minating  laser  output  with  an  acousto-optic  (AO)  device  to  serially  step  through  a 
sequence  of  desired  frequencies.  A  time-frequency  representation  of  a  step-frequency 
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(a) 


detector  element 


Figure  6-1.  Detection  Geometry,  (a)  Filled-aperture  array,  (b)  synthetic-aperture 
array. 


Figure  6-2:  Time-Frequency  Representation  of  Step-Frequency  Illumination  Pulse. 

illumination  is  shown  in  Figure  6-2.  We  require  that  the  duration  of  a  frequency 
step,  At,  be  long  enough  so  that  returns  from  near  and  far  points  in  the  target  can 
simultaneously  have  the  same  frequency. 


At  >  2Zohj/c  ,  (6.5) 

where  Zabj  is  the  depth  of  the  object  and  c  is  the  speed  of  light.  The  duration  of  a 
step-frequency  sequence  will  be 


T  =  NAt  , 


(6.6) 


where  N  is  the  number  of  individual  illuminating  frequencies  or,  equivalently,  the 
number  of  range  bins  desired  in  the  three-dimensional  image. 

There  is  also  an  upper  limit  on  the  speed  of  the  speckle  translating  at  the  detector 
array  that  can  be  tolerated.  We  require  that  a  speckle  move  no  more  than  a  half  of  a 
speckle  diameter  during  an  entire  step-frequency  sequence.  This  requirement  insures 
that  the  angular  sampling  in  the  synthetic-aperture  direction  will  meet  the  Nyquist 
criterion  for  each  of  the  frequencies: 


vs  < 


P 

2  NAt  ' 


(6.7) 
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Substituting  in  expressions  for  vs ,  p,  and  Ai,  we  find  that  this  criterion  restricts  the 
angular  velocity  of  the  target  that  can  be  tolerated: 


Ac 

8  ND0bjZ0bj 


(6.8) 


It  is  worth  noting  that  this  expression  is  independent  of  range.  As  a  sample  calcula¬ 
tion,  consider  an  object  for  which  D0bj  =  5m  and  Z0bj  =  2m.  Suppose  also  that  we 
desire  N  =  10  and  that  we  are  imaging  at  A  =  1  pm.  The  maximum  angular  velocity 
that  can  be  tolerated  for  such  an  object  is  found  to  be  u  <  .375  radians/second. 
This  sample  calculation  suggests  that  the  requirement  on  object  angular  frequency  is 
not  particularly  restrictive  and  that  there  are  likely  to  be  many  objects  that  can  be 
imaged  with  serial-frequency  illumination.  Objects  that  rotate  at  a  rate  faster  than 
the  criterion  prescribed  by  Eq.  (6.8)  will  give  data  that  is  undersampled  (aliased)  in 
one  angular  dimension.  It  is  possible  that  the  opacity  constraint  could  perform  in¬ 
terpolation  in  addition  to  phase  retrieval  to  accommodate  such  undersampled  data. 
Of  course  this  conjecture  would  need  to  be  validated  and  quantified.  Note  that 
restrictions  on  target  angular  velocities  are  relaxed  as  the  array  is  built  up  in  the 
synthetic-aperture  direction.  In  the  case  of  a  completely  filled  array,  only  a  single 
step-frequency  illumination  sequence  is  needed  to  collect  the  entire  data  set  and  the 
three-dimensional  Fourier  volume  will  be  mildly  skewed  in  the  illumination-frequency 
dimension  for  targets  that  rotate  particularly  fast. 


6.5  PARALLEL-FREQUENCY  ILLUMINATION 

In  serial-frequency  illumination,  speckle  measurements  for  differing  frequencies  are 
differentiated  by  temporal  multiplexing.  An  alternative  is  to  illuminate  the  object 
with  all  frequencies  simultaneously  and  use  a  dispersive  element  prior  to  detection  to 
differentiate  frequencies.  This  mode  of  operation  is  referred  to  as  parallel-frequency 
illumination  and  is  depicted  in  Figure  6-3. 

Parallel-frequency  illumination  is  attractive  since  an  AO  device  is  not  required  to 
produce  a  precisely  timed  illuminating  waveform,  as  is  true  in  the  serial- frequency 
illumination  mode.  Multiple  frequencies  can  be  created  by  operating  the  laser  at 
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Figure  6-3.  Parallel-Frequency  Illumination  Mode 
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multiple  longitudinal  modes  simultaneously.  On  the  other  hand,  the  resolving  power 
needed  to  differentiate  frequencies  with  a  dispersive  element  is  demanding.  The  fre¬ 
quency  interval  required  to  meet  the  Nyquist  criterion  is  given  by  [6.1] 

<  rf-.  (6-9) 

£^obj 

The  resolving  power  of  a  dispersive  element  is  defined  as 


_A_ 

AA 

v 

zV  ’ 


(6.10) 

(6.11) 


where  A  and  v  represent  the  center  wavelength  and  frequency,  respectively.  Substi¬ 
tuting  in  the  expression  for  the  maximum  frequency  interval  (Eq.  (6.9)),  we  have 


n  > 
> 


2Z0bjU 

c 

2Z0bj 

A  ' 


(6.12) 

(6.13) 


For  Z0ij  =  2m  and  operating  at  A  =  1  fxm  we  have  bu  <  75 MHz  and  71  >  4  x  106. 
The  resolving  power  of  a  grating  is  given  by  [6.5] 


7 Zg  =  mM  , 


(6.14) 


where  m  is  the  order  and  M  is  the  number  of  grating  lines.  A  blazed  grating  design 

could  be  used  to  direct  the  light  to  a  higher  order - perhaps  the  5th  order.  In 

order  to  achieve  a  resolving  power  of  1Z  =  4  x  106  we  require  M  =  8  x  105  in  a 
diffraction-limited  grating.  With  a  grating-line  separation  of  1  fim,  such  a  grating 
would  be  80  cm  in  length.  Such  an  element  would  be  impractical  and  motivates 
us  to  investigate  other  possibilities.  Some  spectrograph  designs  rely  upon  multiple 
passes  through  a  single  grating,  which  may  make  a  grating  design  more  feasible; 
but  this  would  require  additional  study.  For  a  proof-of-principal  experiment  using 
small  objects  on  a  1  km  test  range,  a  grating  would  be  feasible.  For  example,  if  we 
are  making  a  target  of  depth  2  cm,  then  the  grating  would  have  to  be  just  8mm  in 
length. 
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A  Fabry-Perot  interferometer  can  be  used  to  achieve  resolving  powers  that  are 
significantly  larger  than  gratings.  The  expression  for  the  resolving  power  for  a  Fabry- 
Perot  interferometer  is  given  by  [6.6] 


'R./p 


2n  ndr 
A(1  —  r2)  ’ 


(6.15) 


where  n  is  the  index  of  refraction  of  the  material  between  reflective  surfaces,  d  is  the 
separation  between  reflective  surfaces,  and  r  is  the  reflectivity.  As  an  example,  if 
we  select  r  =  .95,  n  =  1,  and  d  —  ,1m,  the  resolving  power  is  found  to  be  7 Zjp  = 
6.1  x  106,  which  is  sufficient  for  the  previous  example.  Unfortunately,  the  Fabry- 
Perot  interferometer  only  passes  one  narrow  spectral  band  at  a  time  and  rejects  the 
remaining  light.  This  inefficiency  would  not  be  acceptable  for  our  imaging  application 
in  which  photons  are  a  highly  valued  commodity.  One  can  imagine  a  scheme  in  which 
the  rejected  light  is  repeatedly  recycled  through  the  Fabry-Perot  interferometer  at 
slightly  different  angles  on  each  pass  to  yield  a  high-efficiency  system.  Such  a  device 
would  take  further  study. 


The  differentiation  of  frequencies  can  also  be  accomplished  with  heterodyne  de¬ 
tection.  In  this  mode  of  operation  it  is  natural  to  try  to  detect  both  Fourier  modulus 
and  phase  information.  If  the  phase  information  is  corrupted  due  to  instabilities  in 
the  local  oscillator  or  mechanical  instabilities  in  the  optical  distribution  of  the  local- 
oscillator  signal  to  each  detector  in  the  array,  then  a  phasing  algorithm  (see  Section  4) 
will  be  needed.  Opacity  can  also  help  the  phasing  algorithm. 


6.6  CONCLUSIONS 

Various  embodiments  of  three-dimensional  imaging  with  frequency-diverse  speckle 
have  been  described.  The  serial-frequency  illumination  imaging  mode  offers  a  simple 
hardware  design  and  appears  to  be  valid  for  a  significant  range  of  objects  of  interest. 
This  mode  relies  heavily  upon  an  opacity  constraint  to  accomplish  phase  retrieval.  In 
the  parallel-frequency  illumination  imaging  mode,  a  practical  dispersive  element  with 
sufficient  resolving  power  and  high  efficiency  has  not  been  identified  for  the  case  of 
large  objects.  Opacity  can  also  play  a  role  with  phasing  algorithms  when  heterodyne 
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detection  is  employed.  In  any  case,  opacity  promises  to  play  a  role  in  post-detection 

processing  of  the  data  for  phase  retrieval,  bandwidth  extrapolation,  or  interpolation. 

Appendix  C  is  a  reprint  of  a  publication  showing  the  first  successful  reconstruc¬ 
tions  of  a  3-D  image  from  computer-simulated  data  using  the  PROCLAIM  technique. 
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7  LABORATORY  DEMONSTRATION  OF  QUASI-HETERODYNE 

IMAGING 


We  conducted  a  laboratory  experiment  to  demonstrate  quasi-heterodyne  imaging 
and  correction  of  phase  errors.  The  general  approach  was  to  record  pupil  plane 
complex  optical  field  data  by  using  a  detector  array.  We  also  introduced  quadratic 
phase  errors  by  imparting  defocus  to  the  reference  beam.  We  obtained  multiple 
realizations  of  the  pupil  plane  field  by  rotating  the  target  slightly. 


Target  on  Rotation  Stage 


Figure  7-1.  Experimental  Setup  Used  for  Laboratory  Data  Collection. 


7.1  LABORATORY  DATA  COLLECTION 

The  diagram  for  the  experimental  setup  is  shown  in  Fig.  7-1.  Light  from  an  argon 
laser  is  collimated  and  separated  into  object  and  reference  beams.  The  object  beam 
illuminates  a  target  which  consists  of  cutout  pieces  of  retro-reflective  tape  mounted  on 
glass.  The  target  occupies  an  area  that  is  roughly  2  cm  x  2  cm.  Light  reflected  by  the 
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target  passes  back  through  a  beamsplitter  and  then,  without  imaging  optics,  to  a  CCD 
detector  array.  The  reference  beam  of  light  travels  through  the  beam  splitter  and  to  a 
lens  that  creates  a  point  source  that  is  initially  at  the  same  distance  from  the  CCD  array 
as  the  object,  and  is  located  to  coincide  with  the  center  of  the  target.  Light  from  this 
reference  beam  then  propagates  to  a  mirror  that  is  mounted  on  a  piezo-electric 
transducer.  Longitudinal  motion  of  the  transducer  is  used  to  change  the  phase  of  the 
reference  beam.  The  reference  beam  then  propagates  through  a  beamsplitter  to  the 
detector  array  that  records  the  interference  pattern  of  the  object  and  reference  beams. 

The  data  recording  procedure  consists  of  recording  four  interferograms  for  each 
realization  of  the  pupil  plane  field.  These  four  interferograms  are  separated  in  phase  by 
90°,  which  we  accomplish  by  motion  of  the  piezo  mirror.  In  this  manner  we  can  record 
the  complex  values  of  the  optical  field  [7.1].  To  obtain  independent  realizations  of  the 
pupil-plane  field,  we  rotate  the  object  is  rotated  so  that  new  speckles  occupy  the  CCD 
array. 

For  the  initial  set  of  data,  we  recorded  a  series  of  64  independent  realizations  of  the 
pupil-plane  field.  These  frames  of  data  are  initially  stored  on  a  personal  computer  and 
are  then  passed  to  a  network  of  workstation  computers  for  processing. 

7.2  IMAGE  RECONSTRUCTION 

From  a  set  of  four  interferograms  we  compute  the  complex  values  of  the  optical 
field  in  the  aperture  (CCD)  plane.  Then  by  Fourier  transformation  we  compute  a 
complex-valued  image.  Figure  7-2  shows  a  noncoherent  average  of  ten  of  the  image 
realizations  when  no  phase  error  was  introduced. 
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Figure  7-2.  Image  from  Quasi-Heterodyne  Laboratory  Experiment.  The  average  of  the 
intensities  of  ten  reconstructed  images  is  show. 

7.3  SECTION  7  REFERENCES 
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8  DETERMINISTIC  SPARSE  ARRAY  DESIGN 


As  discussed  in  Section  3,  given  a  large  enough  number  of  speckle  realizations,  we 
can  reconstruct  a  diffraction-limited  incoherent  (speckle-free)  image  of  a  laser- 
illuminated  target.  But  to  do  so  from  a  sparse  aperture,  we  must  design  the  pattern  of 
subapertures  so  that  the  autocorrelation  function  of  the  aperture  fills  optical  transfer 
function  (OTF)  space.  If  there  are  gaps  in  OTF  space,  then  we  cannot  recover  those 
missing  spatial  frequencies  by  Wiener  filtering.  At  the  same  time,  to  minimize  cost  we 
may  wish  to  use  the  least  number  of  subapertures  possible.  Hence  we  have  the  problem 
of  finding  a  pattern  with  the  minimum  number  of  subapertures  whose  autocorrelation 
forms  a  filled  pattern.  In  this  section  we  describe  two  such  patterns,  one  of  which,  the 
“arrow,”  fills  a  square  OTF  and  the  second  of  which,  the  “box  +  2,”  fills  a  square  OTF 
with  two-fold  redundancy.  Unlike  other  subaperture  pattern  designs,  they  are 
completely  deterministic  and  are  immediately  scalable  to  any  array  size.  We  wish  to 
define  the  patterns  on  square  array  because  that  matches  the  format  of  the  image 
formation/reconstruction  processing  we  perform  in  the  computer. 

Before  showing  the  more  efficient  arrow  pattern,  we  consider  the  pattern  shown  in 
Figure  8-1,  which  has  the  shape  of  the  letter  “U”  turned  on  its  side.  (Rotating  the 
pattern  by  90°  gives  another  pattern  that  also  fills  OTF  space,  so  the  pattern  that  looks 
like  an  upright  “U”  has  the  same  properties.)  For  a  subaperture  pattern  defined  on  a 
square  grid  with  L  positions  on  a  side,  it  produces  a  filled  autocorrelation  of  size 
(2L  -  1)  x  (2L  -  1)  using  (3L  -  2)  subapertures.  For  the  example  shown  in  Figure  8- 
1,  for  a  subaperture  pattern  with  L  =  8,  it  forms  a  square  autocorrelation  having 
(2L  -  1)  =  15  points  on  a  side  from  (3L  -  2)  =  22  subapertures. 

Figure  8-2  shows  an  example  of  the  arrow  array,  which  is  the  most  efficient  array 
we  have  found  to  date.  The  arrow  array  is  defined  on  a  square  grid,  which  has  the  first 
row  filled,  the  first  column  filled,  and  the  diagonal  between  them  filled,  except  for  the 
(2,  2)  point,  which  is  missing.  (The  array  can  also  be  rotated  by  90°.)  It  yields  a  filled 
autocorrelation  of  size  (2L  -  1)  x  (2L  -  1)  using  only  (3L  -  3)  subapertures,  one  less 
than  the  U  array.  For  the  example  shown  in  Figure  8-2,  it  forms  a  square 
autocorrelation  of  (2L  -  1)  =  15  points  on  a  side  from  only  (3L  -  3)  =  21  subapertures. 
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Figure  8-1.  U  Array  Pattern.  Figure  8-2.  Arrow  Array  Pattern. 

Table  8-1  compares  the  number  of  subapertures  in  the  U  and  arrow  patterns  with 
those  designed  by  a  genetic  algorithm  [8.1].  Of  the  five  cases  given  in  Figure  4  of  [8.1], 
the  U  array  ties  the  genetic  algorithm  in  two  cases  and  beats  it  in  one.  The  arrow  array 
ties  the  genetic  algorithm  in  two  cases  and  beats  it  in  three. 

Table  8-1.  Number  of  Subapertures  to  Fill  (2L-l)x(2L-l)  Autocorrelation  Space. 


Array  Size 
(L) 

Autocorr. 
Size  (2L-1) 

Genetic 

U 

Arrow 

4 

7 

9 

10 

9 

6 

11 

16 

16 

15 

7 

13 

19 

19 

18 

8 

15 

21 

22 

21 

16 

31 

54 

46 

45 

Figure  8-3,  taken  from  Reference  [8.1]  shows  graphically  that  the  Arrow  array 
significantly  outperforms  the  computer-designed  arrays  for  large  array  sizes. 

Depending  on  the  desired  effect,  other  pattern  design  problems  can  be  defined.  For 
example,  what  pattern  with  the  minimum  number  of  subapertures  yields  an  OTF  with 
two-fold  redundancy?  This  would  give  twice  the  signal-to-noise  ratio  as  a  pattern 
simply  filling  OTF  space.  We  designed  the  “box  +  2”  pattern,  an  example  of  which  is 
shown  in  Figure  8-4.  It  yields  a  filled,  doubly  redundant  autocorrelation  of  size 
(2L  -  1)  x  (2L  -  1)  using  only  (4L  -  2)  subapertures.  For  the  example  shown,  L  =  9, 
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(2L  -  1)  =  17,  from  only  (4L  -  2)  =  34  subapertures,  compared  with  36  for  a  pattern 
designed  with  the  genetic  algorithm. 

8.1  SECTION  8  REFERENCES 


[8.1]  P.D.  Henshaw  and  N.R.  Guivens,  Jr.,  “Genetic  Algorithms  for  Unconventional 
Imaging,”  in  Inverse  Optics  III.  Proc.  SPIE  2241.  257-265  (April  1994). 


Figure  8-3.  Comparison  of  Number  of  Subapertures  Needed  to  Synthesize  a  Filled 
OTF  as  a  Function  of  OTF  (Autocorrelation)  Size. (From  [8.1].) 


Figure  8-4.  Box+2  Array  Pattern. 
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9  INTERMEDIATE-PLANE  DETECTION 


Usually  conventional  imaging  sensors  make  measurements  in  a  focal  plane  of  an 
optical  system,  whereas  many  unconventional  imaging  sensors,  including  the  ones 
considered  in  this  report,  make  measurements  in  the  pupil  plane.  An  exception  is  the 
case  in  which  each  subaperture  of  the  receiver  is  small  compared  with  the  diameter  of  a 
pupil-plane  laser  speckle.  If  we  put  an  array  of  detectors  in  such  a  focal  plane,  then  the 
intensity  measured  by  any  one  CCD  element  (or  any  collection  of  them)  is  equivalent  to 
the  intensity  in  the  pupil  plane.  However,  each  CCD  element  then  “sees”  laser  light 
from  only  a  small  portion  of  the  sky.  Consequently,  a  CCD  array  in  the  focal  plane  of  a 
small  subaperture  can  (1)  rapidly  steer  the  field-of-view,  (2)  have  a  large  instantaneous 
field-of-regard,  (3)  gather  data  on  multiple  well-separated  target  complexes  in  parallel, 
and  (4)  suppress  background  light.  These  advantages  of  a  CCD  array  in  the  focal  plane 
of  each  subaperture  must  be  weighed  against  some  disadvantages,  including  the 
relatively  slow  read-out  time  associated  with  most  CCD  arrays  as  opposed  to  an  array 
of,  say,  avalanche  photo-diodes. 

For  those  situations  in  which  a  CCD  array  in  the  back  focal  plane  is  desirable,  we 
can  obtain  additional  design  options  by  a  new  approach  that  entails  putting  the  CCD 
array  in  a  plane  that  is  not  the  focal  plane  of  the  subarray  optics.  Specifically,  it  allows 
the  use  of  subapertures  that  are  larger  in  diameter  than  the  speckle  size  in  the  pupil 
plane;  that  is,  the  target  can  be  partially  resolved  by  the  subaperture.  This  would  allow 
a  system  that  can  image  targets  at  a  wider  range  of  distances. 

This  situation  is  illustrated  in  Figure  9-1.  In  the  figure,  the  first  plane  behind  the 
lens  represents  the  focal  plane  of  the  subaperture.  Parallel  rays  from  the  target  come  to 
focus  in  that  plane.  The  second  plane  behind  the  lens  represents  the  new  proposed 
position  for  the  array  of  detectors.  Rays  striking  a  point  in  that  plane  come  from  the 
same  general  direction  as  the  rays  striking  the  corresponding  (demagnified)  position  in 
the  focal  plane.  However,  the  rays  coming  to  the  second  plane  arrive  from  a  wider 
range  of  angles  than  rays  hitting  the  corresponding  point  in  the  focal  plane.  Hence  there 
is  less  field-of-view  selectivity  in  the  second  plane.  However,  notice  that  the  rays  hitting 
that  point  all  pass  through  a  single  point  in  front  of  the  pupil,  in  a  plane  indicated  in  the 
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figure  by  the  second,  longer  line  a  distance  in  front  of  the  lens.  In  effect,  the  optical 
system  is  forming  an  image  of  the  speckle  pattern  that  is  in  that  plane  in  front  of  the 
lens.  Hence,  if  the  lens  has  a  diameter  of  more  than  one  speckle,  the  detector  can  record 
multiple  speckles.  As  seen  from  the  lower  part  of  the  figure,  rays  hitting  different 
pixels  come  from  different  parts  of  the  lens,  i.e.,  from  different  speckles. 


Figure  9-1.  Detector  Array  in  a  Plane  Behind  the  Focal  Plane  of  the  Lens. 

(a — upper)  Rays  arriving  at  a  pixel  pass  through  a  single  point  in  a  plane  in  front  of  the 
lens;  (b — lower)  rays  hitting  different  pixels  come  from  different  parts  of  the  lens,  i.e., 
from  different  speckles. 


In  summary,  for  a  detector  array  in  focal  plane,  each  pixel  sees  a  different  IFOV, 
thereby  maximizing  spatial  background  suppression,  and  each  pixel  integrates  over  the 
entire  entrance  aperture.  Then  we  must  have  subaperture  diameter  less  than  a  speckle 
size.  At  the  other  extreme,  for  a  detector  array  in  a  plane  that  re-images  the  aperture 
plane,  each  pixel  sees  the  FOV  of  the  telescope,  yielding  minimum  background 
suppression,  and  each  pixel  sees  a  small  portion  of  the  entrance  aperture,  and  so  we  can 
measure  many  speckles  over  a  subaperture.  Both  of  these  extremes  are  poor  for  the 


case  of  a  small  number  of  speckles  per  subaperture.  The  new  method  of  putting  the 
detector  array  in  an  intermediate  plane,  behind  the  focal  plane,  then  each  pixel  sees  a 
somewhat  larger  IFOV  than  for  the  focal  plane,  yielding  partial  background 
suppression,  and  each  area  of  the  detector  sees  a  different  aperture-plane  speckle, 
allowing  multiple  speckles  to  be  measured  simultaneously.  Hence  it  allows  a  trade-off 
of  IFOV  (background  suppression)  against  the  ability  to  measure  multiple  aperture- 
plane  speckles  simultaneously. 
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APPENDIX  A 


This  appendix  describes  the  theory  behind  most  of  the  imaging  modalities  analyzed 
in  this  report.  It  is  from  J.R.  Fienup,  “Unconventional  Imaging  Systems:  Image 
Formation  from  Non-imaged  Laser  Speckle  Patterns,”  in  Emerging  Systems  and 
Technologies.  S.R.  Robinson,  Ed.,  Vol.  8  of  The  Infrared  &  Electro-optical  Systems 
Handbook,  (ERIM,  Ann  Arbor,  1993),  Chapter  1.5,  pp.  92-109. 

A  reader  unfamiliar  with  pupil-plane  laser  imaging  modalities  would  benefit  from 
reading  this  appendix  first. 
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1.5  IMAGE  FORMATION  FROM  NONIMAGED  LASER  SPECKLE  PATTERNS 

Written  by  James  R.  Fienup 

In  most  conventional  imaging  systems,  an  image  is  formed  by  an  optical  system 
and  the  electromagnetic  radiation  (light,  when  at  optical  wavelengths)  is  de¬ 
tected  in  the  image  plane  directly.  However,  if  the  scene  or  object  being  imaged 
is  illuminated  by  coherent  radiation,  as  from  a  laser,  then  several  alternative 
approaches  can  be  used  to  form  an  image  without  the  use  of  image-forming 
optics.  Images  can  be  formed  from  data  collected  in  a  plane  where  the  aperture 
of  an  imaging  system  would  ordinarily  be  located.  This  is  possible  because  the 
electromagnetic  field  in  the  plane  of  the  receiver  is  related  to  the  field  reflected 
by  the  scene  by  an  invertible  transform.  The  transform  is  the  Fresnel  transform 
if  the  scene  is  in  the  near  field  of  the  receiver,  and  the  Fourier  transform  if  in 
the  far  field.393  If  field  measurements  are  made  at  an  array  of  points  in  an 
aperture  plane  (i.e.,  without  a  focusing  system)  by  heterodyne  or  holographic 
methods,  then  the  image  formation  can  be  accomplished  by  an  inverse  Fourier 
(or  Fresnel)  transform  in  a  digital  computer.  This  is  commonly  done  for  mi¬ 
crowave  synthetic  aperture  radar  (SAR)  systems. 

The  symbols  used  in  this  section  are  defined  in  Table  1.2. 

The  resolution  limit  imposed  by  diffraction  is 


where  p  =  resolution  at  the  object,  A  —  wavelength,  R  —  range  to  target,  and 
Da  =  diameter  of  the  aperture  or  collecting  array. 

The  familiar  1.22 \R/Da  is  replaced  by  \R/Da  because  we  are  considering  a 
square  aperture  of  length  Da.  To  image  distant  objects  with  very  fine  resolution, 
this  limit  on  resolution  forces  the  use  of  short  wavelengths  in  order  to  allow 
for  apertures  of  practical  size.  However,  for  shorter  wavelengths,  particularly 
in  the  ultraviolet  through  the  midinfrared,  heterodyne  measurements  of  the 
field  are  very  difficult.  Furthermore,  at  these  wavelengths  the  turbulence  of 
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Table  1.2  Symbols,  Nomenclature,  and  Units 


Symbol 

Nomenclature 

Units 

c 

Speed  of  light 

m/s 

C 

Autocovariance  of  aperture-plane  intensity 

W2/m2 

da 

Detector  spacing 

m 

daT\d 

Area  of  one  detector 

m2 

Da 

Diameter  of  receiver  aperture 

m 

Do 

Object  diameter 

m 

f 

Optical  field  in  object  plane 

VW/m2 

ft 

Incoherent  object  brightness 

W/m2 

F 

Optical  field  in  aperture  plane 

VW/m2 

Fi 

Fourier  transform  of  incoherent  object 

W/m2 

g 

Optical  field  in  image  plane 

VW/m2 

h 

Planck’s  constant 

J  s 

i 

v~ 

Lc 

Coherence  length  of  laser 

m 

Lp 

Length  of  laser  pulse 

m 

M 

Number  of  detectors  per  aperture-plane  speckle 

flbp 

Number  of  resolved  bright  image  points  at  a  given  resolu¬ 
tion 

nd 

Average  number  of  photons  per  detector 

ns 

Number  of  speckles  in  the  aperture 

.  N 

Number  of  independent  realizations  of  the  speckle  pattern 

Nm 

Maximum  number  of  independent  realizations  available 

P 

Coherent  impulse  response 

P 

Pupil  function 

r 

Autocorrelation  of  coherent  object 

W/m 

n 

Autocorrelation  of  incoherent  object 

W2/m2 

ra 

Object  average  reflectivity 

R 

Range  to  object 

m 

s 

Point  spread  function,  |p|2 

s 

(Normalized)  optical  transfer  function 

So 

(Unnormalized)  optical  transfer  function 

SBP  max 

Maximum  1-D  space-bandwidth  product 

SNR 

Signal-to-noise  ratio 

( u,v ) 

Spatial-frequency  coordinates 

m 

Ur 

(u2  +  v2)n 

m 

Va 

Speckle  velocity  in  aperture  plane 

m/s 

Vo 

Object  velocity 

m/s 

(x,y) 

Image-plane  coordinates 

m 

(Xaja) 

Aperture-plane  coordinates 

m 

(Xojo) 

Object-plane  coordinates 

m 

Greek: 

Aii 

Time  for  pattern  to  move  across  one-quarter  of  the  speckle 
width 

s 
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Table  1.2  (continued) 


Symbol 

Nomenclature 

Units 

Time  for  pattern  to  move  across  aperture 

s 

Time  before  object  rotates  too  much 

s 

A z 

Object  depth  along  line  of  sight 

m 

T)d 

Detector  duty  cycle 

Quantum  efficiency  of  the  detector 

m 

X 

Wavelength 

P 

Complex  coherence  factor 

m 

Pa 

Speckle  width  at  aperture 

Po 

Potential  resolution  at  object 

m 

Tatm 

One-way  atmospheric  transmittance 

Transmittance  of  receiver  optics  (if  any) 

Fraction  of  light  reflected  in  desired  polarity  - 

rad 

Estimate  of  aperture-plane  phase 

Aperture-plane  phase 

rad 

(D 

Object  rotation  rate  (axis  perpendicular  to  line  of  sight) 

rad/s 

Operators: 

9 

Fourier  transformation 

* 

Convolution 

Autocorrelation 

the  atmosphere  causes  phase  errors  in  the  field  propagated  through  it,  which 
can  severely  limit  the  resolution  of  a  heterodyne  system.  If  it  is  possible 
reconstruct  an  image  from  the  intensity  of  the  field  m  the  aperture  plane  then 
these  limits  to  resolution  caused  by  aberrations  can  be  overcome.  Figure  1.40 
shows  examples  of  configurations  that  can  be  used  to  gather  the  aPe^ure'P^^ 
intensity.  These  laser  intensity  patterns  are  speckled  if  the  scene  is  rough  in 
comparison  to  the  wavelength  of  illumination.  In  this  section,  some  methods 
for  forming  images  from  nonimaged,  aperture-plane  speckle  patterns  are  dis¬ 
cussed  First  some  terms  are  defined  and  a  few  methods  are  briefly  desc 
Then  one  particular  method,  imaging  laser  correlography,  is  described  in  more 
detail,  to  serve  as  an  example  of  a  nonconventional  imaging  method  using 

aperture-plane  speckled  intensity  patterns.  ...  .  ,  ,  f 

These  techniques  generally  require  that  the  scene  be  of  limited  exteut- 
that  it  be  a  bright  object  on  a  dark  background.  They  have  been  developed 
primarily  for  the  applications  of  fine-resolution  imaging  of  space  objects  by 
Lund-based  systems  and  of  discriminating  decoys  from  reentry  vehicles  for 
the  Space  Defense  Initiative  (SDI)  midcourse  phase.  Another  application  could 
be  the  imaging  of  a  ground  scene  from  an  overhead  platform  if  only  a  well- 
defined  portion  of  theground  were  illuminated  by  a  coherent  source.  Relying 
on  powerful  sources  of  short-wavelength,  coherent  radiation,  these  techniques 
have  been  in  development  only  since  the  1960s.  The  first  concerted  effort 
to  develop  these  ideas  began  at  the  Woods  Hole  Summer  Study  in  19  . 
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(a) 


Alternative  detector  array  geometries 


Fig.  1.40  Imaging  geometries:  (a)  The  intensity  of  coherent  radiation  scattered  from  an 
object  is  measured  by  an  array  of  detectors  without  intervening  optics,  (b)  Alternatively,  a 
collection  lens  can  be  used  to  gather  the  light  onto  each  of  the  individual  detectors,  which 
can  have  smaller  areas,  (c)  Another  alternative  is  to  reimage  the  pupil  plane  of  a  telescope 
onto  a  detector  array  of  small  dimensions,  such  as  a  CCD  array. 


1.5.1  Basic  Imaging  Equations 

Because  of  the  mixture  of  coherent  and  incoherent  imaging  concepts  inherent 
in  imaging  correlography,  it  is  necessary  to  define  here  the  basic  imaging 
equations  that  will  be  needed  later  in  this  section. 

The  object  is  assumed  to  be  illuminated  by  radiation  that  is  both  spatially 
and  temporally  coherent  over  the  volume  of  the  object.  Here  we  refer  to  the 
illuminated  scene  as  the  object.  Let  the  complex  field  reflected  by  the  object 
be  f'(x0,y0 )  in  a  plane  perpendicular  to  the  line  of  sight  (LOS).  Then  the  field 
in  the  aperture  plane  at  the  detector  array,  also  for  simplicity  assumed  to  be 
perpendicular  to  the  LOS,  is  given  by  the  Fresnel  transform393: 


1  (  i27t\ 

=  jexp 

177  o  9 . 

XR{X°  + 

X  fJf'(x0,yo)  exp 

ITT  o  2\ 

ysix-  +y°\ 

x  exp 

-  i2-n\ 

XS  U*I» 

+  yay0 ) 

dx0 

(1.145) 
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In  this  section,  all  integrals  are  over  ( -  *,*).  For  the  case  of  a  diffusely  scat¬ 
tering  object,  the  intensity  of  this  aperture-plane  optical  field  is  a  speckle 
pattern.  We  are  assuming  that  in  front  of  the  detector  array  is  a  polarizer  that 
selects  a  single  polarization,  which  is  necessary  for  the  methods  described  here 
to  work.  Hence,  throughout  this  section  we  are  considering  radiation  of  a  single 
polarization,  enabling  us  to  employ  scalar  diffraction  theory. 

If  we  absorb  the  constant  phase  term  -2 the  quadratic  phase  term 
( -  tt /\R  )(x2  +  y2a ),  and  the  factor  i  into  the  phase  of F'(xa,ya)  to  form  F(xa,ya), 
and  absorb  the  quadratic  phase  term  (n/\R  )(x20  +  y2)  into  the  phase  off'(x0,y0) 
to  form  f(x0,y0),  then  we  have  the  simpler  Fourier  transform  relationship 


F(xa,ya )  =  |F(wa)|  exp[i4»(*a,ya)]  =  W(*o,y<>)] 


XR 


ff 

—  i2iT 

1 J  f(xo,y0 )  exp 

jjj  ixax0  +  yay0) 

dx0  dy0 


(1.146) 


between  the  field  f(x0,y0 )  at  the  object  plane  and  the  field  F(xa,ya )  at  the 

OQO 

receiver  plane.  This  Fourier  relationship  is  valid  for  the  near-field  conditions 
normally  allowed  only  for  the  Fresnel  transform.  If  the  complex-valued 
aperture-plane  field,  F(xa,ya),  could  be  measured  without  any  phase  errors, 
then  an  inverse  Fourier  transform  performed  in  a  computer  would  yield  a 
diffraction-limited  complex-valued,  coherent  image.  This  image  would  be  given 
by  the  convolution  (with  appropriate  scaling  of  the  coordinates  to  account  for 
the  image  magnification) 


g(x,y)  =  f*p(x,y) 


j_r 

KRJ. 


f(x0,y0)p(x  -  x0,  y  -  ya)  dx0  dy0 


(1.147) 


where  *  denotes  convolution,  p(x,y)  =  &  _1[F(x0,ya)]  is  the  system  (coherent) 
impulse  response,  and  P(xa,ya )  is  the  pupil  function  (unity  over  the  extent  of 
the  aperture  and  zero  outside).  The  corresponding  point  spread  function  for 
incoherent  imaging  through  the  same  aperture  would  be 

s(x,y)  =  |p(x,y)|2  =  y^ISodi,*;)]  ,  (1.148) 

where  S0(u,v)  is  the  unnormalized  optical  transfer  function,  given  by  the  au¬ 
tocorrelation  of  the  pupil  function: 


So(u,v )  =  P®P(u,v)  =  — 


P(xa,ya)P*(xa 


u,  ya  -  v)  dxa  dya  , 

(1.149) 


and  P*  is  the  complex  conjugate  ofP,  and  (u,v)  is  a  spatial-frequency  coordi¬ 
nate.  The  (normalized)  optical  transfer  function  is  given  by 


S(u,v) 


So(u,V ) 

Wi  ■ 


(1.150) 
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Like  any  coherent  image,  g{x,y)  would  be  speckled,  assuming  that  the  surface 
of  the  object  is  rough  compared  to  the  wavelength  of  the  radiation.397  As  the 
object  translates  or  rotates,  a  new  part  of  the  back -scattered  field  falls  on  the 
receiver;  effectively,  the  object  takes  on  a  different  realization  of  the  phase  in 
the  object  plane.  When  we  wish  to  distinguish  these  different  fields  coming 
from  the  object,  we  denote  the  n’th  realization  of  the  object  field  as  fn(x0,y0 ) 
and  the  corresponding  field  at  the  aperture  plane  as  Fn(xa,ya )•  The  image  of 
each  fn(xo,y0 )  would  have  a  different  realization  of  the  speckle  pattern  across 
it.  If  the  field  in  the  aperture  plane  translates  by  the  width  of  the  aperture 
between  the  time  the  two  images  are  collected,  then,  for  the  case  of  a  diffusely 
scattering  object,  the  two  field  distributions  are  independent  and  the  complex 
images  are  uncorrelated.  Conditions  for  achieving  this  are  discussed  in  Sec. 
I.5.4.I.  An  average  over  many  such  object  realizations  would  give 

1  N 

, fi(x0,y0 )  =  <|/n(*o,yo)|V  =  lim  -  2  I  fn(x0,y0)\2  ,  '  (1.151) 

where  <  •  )„  denotes  the  ensemble  average  over  the  speckle  realizations  or  aspect 
angles  and  N  is  the  number  of  independent  realizations.  This  noncoherent 
average  of  many  coherent  objects  yields  the  ordinary  incoherent  reflectivity 
(or  brightness  distribution)  of  the  object,  which  is  a  real-valued,  non-negative 
distribution.  We  denote  the  Fourier  transform  of  the  incoherent  object,  fi(x0,y0), 
as  Fi(u,v).  The  coherent  images  p*fn{x,y)  are  speckled  versions  of  s*fiix,y ),  the 
image  of  the  incoherent  object.  The  normalized  Fourier  transform. 


M-(u,f) 


Fj(u,v) 
Fi( 0,0)  ’ 


(1.152) 


is  the  complex  coherence  factor  (or  the  complex  visibility  function),398  and 
|m-(u,i;)|2  is  the  normalized  spatial  power  spectrum  (i.e.,  the  normalized  squared 
modulus  of  the  Fourier  transform)  of  the  incoherent  object. 


1.5.2  Imaging  Methods 

The  quantity  that  we  assume  to  be  measured,  and  from  which  we  wish  to 
reconstruct  an  image,  is  the  aperture-plane  intensity,  |F(xa,ya)|2  -  |F'(xa,ya)|2- 
The  square  root  of  the  intensity,  |F(xa,ya)|,  is  the  magnitude  (commonly  called 
the  modulus  or  amplitude )  of  the  field  and  iMxa,ya)  is  the  phase.  Without  the 
phase  of  F{xa,ya),  the  image  cannot  be  computed  simply  by  inverse  Fourier 
transformation.  Since  1 9[f*( -x,-y)]|  =  IWUoOll  =  |F(xa,ya)|,  from  the 
aperture-plane  intensity  it  is  impossible  to  distinguish  between  the  object 
f(x,y)  and  its  twin  f*(-x,-y).  This  twofold  ambiguity  is  usually  acceptable 
since  the  twin  looks  just  like  the  object,  but  rotated  by  180  deg,  and  it  usually 
suffices  for  identifying  the  object.  In  addition,  \&[f(x  -  x\,y  -  yi)]|  =  \&[f{x,y)]\, 
and  so  the  absolute  position  of  the  object  cannot  be  determined  from  the  in¬ 
tensity  |F(xa,ya)|2  alone. 
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We  next  briefly  describe  some  of  the  methods  that  can  be  used  to  form  images 
from  this  seemingly  insufficient  aperture-plane  intensity  data. 

1 .5.2. 1  Holography.  By  inverse  Fourier  transforming  |F(xa,ya)|2,  we  arrive 
at  the  autocorrelation  of  the  object: 

r(x,y)  =  f®ftx,y)  =  - ^jjf(x',y’)f*i.x  -  x',y- y'ldxdy 

=  5'’l|F(x„,yo)|2]  •  <1153) 

If  the  object  consists  of  an  extended  part  of  the  object  plus  a  glint,  separated 
from  the  extended  part  by  at  least  the  width  of  the  extended  part,  then  |i',(xa,ya)| 
is  the  equivalent  of  a  hologram399  of  the  extended  part.  When  this  holographic 
condition  is  met,  the  autocorrelation  consists  of  four  terms,  one  of  which  is  an 
image  of  the  extended  part  of  the  object.  Since  |F(xa,ya)|  is  measured  only 
over  the  aperture  P{xa,ya )  of  width  Da,  the  autocorrelation  and  the  recon¬ 
structed  image  are  limited  in  resolution  (in  object  space)  to  p  =  \R/Da.  This 
imaging  modality  was  analyzed  for  the  case  of  imaging  space  objects  from  the 
earth  by  Refs.  394  and  395. 

1 .5.2.2  Object  Support  from  Autocorrelation  Support.  The  diameter  of  the 
object  f(x,y)  along  any  direction  is  just  half  the  diameter  of  the  autocorrelation 
r(x,y )  along  the  same  direction.  The  support  of  the  object  is  defined  as  a  set 
of  points  (or  an  area)  outside  of  which  the  object  has  a  value  of  zero,  consistent 
with  the  assumption  that  we  are  considering  bright  objects  on  a  dark  back¬ 
ground.  The  support  of  the  autocorrelation  allows  us  to  determine  bounds  on 
the  support  of  the  object.  For  the  sake  of  brevity  these  methods  are  not  detailed 
here,  except  to  say  that  they  involve  taking  the  intersection  of  two  or  more 
translated  versions  of  the  support  of  the  autocorrelation.  The  reader  is  referred 
to  Refs.  401  and  402  for  the  details.  For  favorable  cases,  such  as  a  highly 
nonconvex  object  or  an  object  consisting  of  two  or  more  separated  parts,  the 
upper  bound  on  the  support  of  the  object  can  serve  as  an  approximate  outline 
of  the  object.  This  outline  may  be  sufficient  information  about  the  object  for 
some  applications.  For  less  favorable  cases,  such  as  for  a  convex  object,  the 
upper  bound  only  gives  the  general  size  of  the  support  of  the  object.  In  all 
cases,  the  upper  bound  on  the  support  can  be  used  as  a  support  constraint  for 
phase  retrieval,  as  is  described  in  the  next  section. 

1.5.2. 3  Phase  Retrieval.  An  approach  to  image  formation  from  coherent 
aperture-plane  intensity  is  phase  retrieval.  The  approach  is  to  use  a  computer 
algorithm  to  determine  a  phase  4>(*a,ya),  which,  when  combined  with  |F(xa,ya)| 
from  the  measured  intensity  and  Fourier  transformed,  yields  an  image 

g{x,y)  =  S'  -1{|F(xa,ya)|  exp[i4>(xa,ya)]}  (1.154) 

that  is  consistent  with  the  support  constraint  (i.e.,  has  a  value  close  to  zero  in 
the  region  outside  the  upper  bound  on  the  support  of  the  object).  Dozens  of 
phase  retrieval  algorithms  have  been  proposed— see,  for  example,  Ref.  403. 
Phase  retrieval  algorithms  have  been  particularly  successful  for  real-valued, 
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non-negative  images  such  as  those  one  would  have  for  incoherent  illumination. 
However,  from  \F(xa,ya)\  the  quantity  being  reconstructed  is  the  field,  which 
is  complex  valued.  Then  the  only  constraint  available  is  the  support  constraint. 
In  this  case  it  is  possible  to  reconstruct  an  image,  particularly  for  objects 
consisting  of  two  or  more  separated  parts;  but  it  is  difficult  and  is  not  yet 
reliable  for  general  objects.404,405 

1 .5.2.4  Additional  Measurements.  If  additional  measurements  are  made, 
then  the  problem  of  reconstructing  a  complex- valued  image  becomes  much 
easier.  One  example  for  a  space-based  system  is  to  have  a  diffraction-limited 
aperture  of  moderate  size  embedded  within  a  larger  array  of  detectors.  By  use 
of  a  beamsplitter,  the  intensity  passing  through  the  moderate-sized  aperture 
is  measured  simultaneously  in  two  different  planes:  in  the  focal  plane,  which 
contains  a  diffraction-limited  but  moderate-resolution  image,  and  in  the  plane 
of  the  aperture.406  Then  the  reconstruction  of  the  phase  over  the  diffraction- 
limited  aperture  can  be  accomplished  by  an  iterative  algorithm  that  is  ro¬ 
bust.407,408  This  phase  can  then  be  used  to  help  to  reconstruct  the  phase  over 
the  larger  array,  where  only  intensity  measurements  were  made,  in  a  boot¬ 
strapping  version  of  a  phase  retrieval  algorithm.406 

A  second  approach  using  additional  measurements  is  to  combine  aperture- 
plane  intensity  methods  with  wavefront  sensor  data.  The  imaging  wavefront 
(Hartmann-Idell)  sensor409  combines  wavefront  slope  data  from  a  Hartmann 
array  with  the  aperture-plane  intensity  data  to  reconstruct  the  optical  field, 
from  which  an  image  is  reconstructed  by  Fourier  transformation.  By  tempo¬ 
rally  averaging  the  wavefront  slopes  over  many  speckle  realizations  during  a 
time  interval  less  than  the  atmospheric  correlation  time,  it  is  also  possible  to 
estimate  the  aberrations  due  to  the  turbulent  atmosphere,  which  could  then 
be  subtracted  to  form  a  diffraction-limited  image. 

Other  schemes  using  additional  measurements  are  also  possible. 


1.5.3  Imaging  Correlography 

The  diffraction-limited  autocorrelation,  p*rn(x,y),  which  can  be  computed  from 
the  aperture-plane  intensity  P(xa,ya)\Fn(xa,ya)\2,  is  a  complex-valued,  speckled 
autocorrelation.  As  the  object  translates  or  rotates,  the  aperture-plane  inten¬ 
sity  pattern  moves.  After  sufficient  motion  at  the  receiver  a  new  realization 
of  the  speckled  intensity  pattern  appears,  which  yields  a  new  realization  of 
the  speckled  autocorrelation.  The  diffraction-limited  autocorrelation  of  the 
incoherent  image,  |p|2*r/(:c,y),  can  be  determined  by  averaging  the  squared 
magnitudes  of  the  speckled  autocorrelations  over  many  realizations  of  the 
speckle  pattern410-412: 

1  N 

lim  —  ^  \p*rn(x,y)\2=  6i|p(x,y)]2  +  b2\p\2*ri(x,y)  ,  (1.155) 

N->  x  An  =  i 


where 

b  i  =  &2 


I2 

jf  I fi(x',y'f  dx'  dy' 


(1.156) 
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is  the  square  of  the  total  measured  intensity,  62  is  a  constant  and  N  is  the 
number  of  independent  realizations.  This  result  assumes  a  rough  object  whose 
coherent  autocorrelation  function  is  of  the  form  of  a  delta  function.  The  first 
term  on  the  right-hand  side  of  Eq.  (1.155)  is  a  point-spread  function  of  known 
shape  located  at  the  origin  (the  dc  term)  that  can  be  subtracted,  leaving \p\  *n(xy), 
the  diffraction-limited  autocorrelation  of  the  incoherent  object.  That  is,  the 
noncoherent  average  of  the  autocorrelations  of  the  coherent  images  yields  the 
autocorrelation  of  the  incoherent  image.  Alternatively,  processing  the  same 
intensity  measurements  P{xa,ya)\Fn(xa,yaf  in  a  different  way,  we  compute 
the  autocovariance  of  the  intensity ,  which  is  given  by 


C(u,V )  =  lim  T;  S  [  \p(xa  +  +  v^xa^a) 

N-*  00  Ar  n  =  1  ■' 

X  [| Fn(Xa  +  U,ya  +  v)\2\Fn(xa,ya)\2-  72]  <ha  . 


(1.157) 


where  7  is  the  average  intensity  over  the  aperture  plane.  Using  the  moment 
factoring  theorem  for  circular-complex  Gaussian  fields  (Ref.  414,  p.  44),  it  can 

be  shown  that 


lim  S  [\Fn(xa  +  u,ya  +  v)\2\Fn(xa,yaf  -  I2]  ~  \Fi{u,v)\' 

yv„=i 


N 


(1.158) 


and  so 


C(u,v)  =  SU,tf)|F/(u,u)|2 


(1.159) 


where  S(u,v )  is  the  optical  transfer  function  for  the  aperture.  That  is,  the 
autocovariance  of  the  intensity  of  the  coherent  aperture-plane  field  is  propor¬ 
tional  to  the  squared  magnitude  of  the  Fourier  transform  of  the  incoherent 
object.  This  is  related  to  the  incoherent  autocorrelation  n(x,y)  by  a  Fourier 
transform.  Computation  of  the  autocovariance  of  the  aperture-plane  intensity 
is  completely  equivalent  to  intensity  interferometry,  but  using  realizations 
of  reflected  laser  speckle  patterns  rather  than  incoherent  light  emitted  by  the 

°^Each  pair  of  detectors  over  the  aperture  contributes  to  a  measure  of  C(u,u) 
and  therefore  to  | F/(u,o)|2.  If  the  pair  of  detectors  is  at  locations  (xa\,yai)  an 
(xa2,ya2)  in  the  aperture  plane,  then  they  contribute  to  a  measurement  of 

\Fi(u,v)\2  at  spatial  frequency  (u,v)  -  ±Uci  xa2,  yai  yo2>-  ®  . 

to  the  fact  that  | F^-u^vf  =  \F!(u,v)\2,  and  both  are  measured  simulta¬ 
neously.]  In  radio  interferometry,  the  separations  (u,v)  are  referred  to  as  base¬ 
line  separations.  Redundant  measurements— different  pairs  of  detectors  with 
the  same  baseline  separations— contribute  to  increasing  S(u,v)  and  increasing 

the  SNR  of  the  estimate  of  \Fi(u,v)\  .  .  ,  , 

Having  computed  \Fi(u,v)\  as  described,  an  incoherent  image  can  be  fairly 
reliably  reconstructed  using  a  phase  retrieval  algorithm  even  though  the 
reconstruction  of  the  individual  coherent  images  from  the  intensities  1 1  „  (xa  ,ya  )| 
is  difficult.  The  reconstruction  is  easier  for  incoherent  images  because  m 
addition  to  the  support  constraint,  we  have  the  additional  constraint  that  the 
image  is  a  real-valued,  non-negative  function. 
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Prior  to  image  reconstruction,  it  is  usually  advantageous  to  perform  a  fil¬ 
tering  or  a  weighting  of  \Fi(u,v)\.  This  is  needed  because  the  SNR  for  the 
higher  spatial  frequencies  is  much  poorer  than  that  at  lower  spatial  frequen¬ 
cies  For  the  highest  spatial  frequencies,  the  noise  usually  dominates.  The 
spatial  frequencies  at  which  SNR  <  1  should  be  discarded  or  at  least  greatly 
deemphasized.  This  can  be  accomplished,  for  example,  by  multiplying  |F/(u,t>)| 
by  a  filter  function 

S(u,v)\Fi(u,v)\4  fi  ■jftfi’i 

where  \Fi(u,vf  is  an  estimate  of  \F,(u,v)\2  (which  could  be  an  average  power 
spectrum  for  objects  of  an  appropriate  class)  and  4>(u,u)  is  the  variance  of  the 
noise  in  the  estimate  of  the  autocovariance  function.  This  filter  function  is 
related  to  the  Wiener-Helstrom  filter.416  Other  variations  on  the  filter  function 
have  also  been  used— anything  that  suppresses  the  highest  spatial  frequencies, 
where  SNR  <  1,  can  be  beneficial. 

Figures  1.41  and  1.42  show  an  example  of  a  computer  simulation  of  imaging 
laser  correlography.41^  In  Fig.  1.41  we  see  that  as  more  realizations  of  the 


Fig  1  41  Noncoherent  averages  of  coherent  image  autocorrelations.412  First  row,  left  to 
right:  noncoherent  average  ofN  =  4  coherent  lmage  autocorrelations,  dc  term  and  averaged 
autocorrelation  with  dc  term  subtracted;  second  row:  N  -  32;  third  row:  N  128;  and  fourth 
row:  N  =  1024. 
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Fie  1  42  Image  recovery  from  simulated  correlography  data.41"  (A)  Estimate  of  the  au¬ 
tocorrelation  of  the  incoherent  image  from  .V  =  10.000  coherent  image  autocorrelations; 
(B)  estimate  of  the  Fourier  modulus  of  the  incoherent  object;  iC)  image  reconstructed  1 Yom 
(B)  bv  means  of  the  iterative-transform  phase  retrieval  algorithm;  (D)  Wiener  filter.  (Ei  filtered 
Fourier  modulus;  (Ft  image  reconstructed  from  (El;  <G»  original  incoherent  object;  (HI  Wiener- 
filtered  incoherent  object;  and  ( 1  >  result  of  W  lener  filtering  of  i  C  i. 


speckled  autocorrelations  are  averaged  together,  it  approaches  the  autocor¬ 
relation  of  the  incoherent  object.  Figure  1.42  shows  an  example  of  image  re¬ 
construction  from  simulated  correlography  data  using  the  iterative  transform 
phase  retrieval  algorithm.408  Filtering  the  data  is  necessary  to  suppress  noise 
at  the  higher  spatial  frequencies.  Figure  1.42(F),  which  shows  the  image  re¬ 
constructed  from  the  filtered  data,  compares  favorably  with  Fig.  1.42(H),  which 
shows  the  ideal  image  passed  through  the  same  filter.  Figure  1.43  shows  an 
example  of  imaging  laser  correlography  in  a  laboratory  experiment. 
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Fig.  1.43  Example  of  reconstruction  from  data  gathered  in  the  laboratory.  Top  left:  one 
realization  of  the  measured  speckle  intensity  patterns:  top  right:  incoherent  image  of  the 
object:  bottom  left:  autocorrelation  estimate  computed  using  100  realizations  of  the  speckle 
pattern;  and  bottom  right:  the  reconstructed  image.  <  Image  courtesy  of  J.  D.  Gonglewski, 
P.  S.  Idell,  and  D.  G.  Voelz.  Phillips  Laboratory! 


1.5.4  System  Parameters 

With  the  basic  method  of  imaging  laser  correlography  in  hand,  we  now  describe 
the  requirements  on  the  system— the  transmitter  (laser)  and  the  receiver  (de¬ 
tector  array  ).  For  the  sake  of  an  example  calculation,  we  use  the  parameters 
given  in  Table  1.3. 

1.5.4. 1  Receiver  Requirements.  To  obtain  a  resolution  at  the  object  of  p0 
( =  0.2  m  for  the  example),  an  aperture  with  a  diameter  of  at  least 


Da  = 


KR 

Po 


(1.161) 
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Table  1.3  Example  of  Parameters  for  an  Imaging  Scenario 


Symbol 

Example  Value 

Description 

Do 

5  m 

Object  diameter 

Az 

5  m 

Object  depth  along  LOS 

0) 

1  rad/s 

Object  rotation  rate  (axis  perpendicular  to  LOS) 

Vo 

5  km/s 

Object  velocity  perpendicular  to  LOS 

X 

1  pm 

Wavelength 

R 

1000  km 

Range  to  object 

Da 

5  m 

Receiver  aperture  diameter 

Po 

0.2  m 

Resolution  at  object  =  kRIDa 

Po 

0.2  m 

Speckle  size  in  aperture  =  \R/D0 

( =  5  m  for  the  example)  is  required.  Limits  on  resolution  that  have  to  do  with 
the  SNR  will  be  discussed  later. 

A  fundamental  quantity  needed  to  explain  the  system  requirements  is  the 
characteristic  diameter  pQ  of  the  speckles  in  the  detected  intensity  pattern, 
\Fni.xa,ya)\2.  It  is  given  by 


(1.162) 


( =  0.2  m  for  the  example).  The  maximum  possible  space-bandwidth  product 
of  the  image  in  each  dimension  is 


SBP 


max 


D0  D0Da  _  Da 

Po  Pa 


(1.163) 


( =  25  for  the  example).  This  is  the  maximum  number  of  resolution  elements 
across  the  object,  which  is  equal  to  the  number  of  speckles  across  the  aperture. 
The  actual  number  of  resolution  elements  across  the  object  may  be  fewer 
because  the  resolution  can  be  limited  by  noise,  as  is  described  later.  In  order 
for  |F„(xa,ya)|2  to  be  measured  in  the  aperture  plane  without  aliasing  (to  avoid 
undersampling),  the  center-to-center  spacing  of  the  detectors  must  satisfy 


£a  =  (1.164) 

a  “  2  Wo 

(=  0.1  m  for  the  example),  which  gives  two  samples  per  speckle  in  each  di¬ 
mension.  For  a  filled  array,  DJda  ( =  50  for  the  example)  detectors  are  needed 
in  each  dimension.  However,  since  the  data  used  for  reconstruction,  \Fiiu,v)\  , 
are  given  by  the  autocorrelation  of  the  measured  data,  |F n(xa>ya)|  >  sampling 
of  the  measured  data  on  a  regular  filled  grid  is  not  necessary.  In  fact,  sparse 
arrays  of  detectors  may  be  used.413  The  effective  transfer  function  of  the  system 
is  the  autocorrelation  of  the  aperture  function,  which  is  equivalent  to  the 
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optical  transfer  function.  Issues  regarding  filling  of  the  (u,v)  plane  (Fourier 
aperture)  are  similar  to  those  encountered  in  radio  interferometry.4 

The  intensity  pattern  moves  across  the  aperture  plane,  due  to  the  object’s 
rotation,  at  a  velocity  2 wi?  ( =  2000  km/s  for  the  example),  and  so  it  moves  by 
a  speckle  diameter  in  time  pa/(2ooi? )  =  \/(2u) D0)  ( =  100  ns  for  the  example). 
Therefore  the  integration  time  for  a  single  speckle  pattern  (without  some  form 
of  speckle  tracking)  should  be  limited  to  about 
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( =  25  ns  for  the  example).  The  intensity  pattern  also  moves  across  the  aperture 
plane  at  a  velocity  2v0  (=  10  km/s  for  the  example)  due  to  the  velocity  of  the 
object.  Note  that  at  long  ranges  the  velocity  of  the  aperture-plane  speckles  due 
to  object  rotation  are  much  greater  than  that  due  to  object  velocity,  except  for 
very  slow  rotation  rates.  If  the  detectors  have  an  integration  time  or  a  readout 
time  that  is  long  compared  with  A<i,  then  the  detectors  must  be  gated  over 
an  interval  of  time  Mi  (or  the  speckles  must  be  tracked).  For  a  pulsed  laser 
of  pulse  length  Lp ,  a  long-readout  detector  could  only  use  the  fraction  cM\!Lp 
of  the  energy  of  the  pulse,  where  c  is  the  speed  of  light. 

The  time  it  takes  the  speckles  to  pass  over  the  receiver  aperture  of  diameter 
Da  is 


Ma 


Da 

2  ojR 


(1.166) 


( =  2.5  |xs  for  the  example).  This  is  the  time  necessary  to  obtain  an  independent 
realization  of  the  speckle  pattern  and  defines  the  maximum  rate  at  which 
independent  realizations  can  be  collected  (400,000  independent  realizations 
per  second,  for  the  example).  The  maximum  integration  time  for  a  single 
speckle  realization  is  less  than  Ma  even  with  speckle  tracking.  Since  there 
are  SBPmax  =  D0I po  ( =  25  for  the  example)  diffraction-limited  resolution 
elements  across  the  object,  a  rotation  of  the  object  by  roughly  A RI(D0Da)  ra¬ 
dians  will  cause  the  object  brightness  distribution  to  begin  to  change  sub¬ 
stantially.  This  would  occur  in  time 
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( =  40  ms  for  the  example).  Unless  multiple  revolutions  of  the  object  are  used, 
this  is  also  the  maximum  amount  of  time  available  for  gathering  an  image  of 
the  object  for  a  single  aspect  angle.  During  time  Mm,  the  maximum  number 
of  independent  realizations  available  is 


Mm  2 \R2 
Ma  D0Dl 


(1.168) 


(  =  16,000  realizations  for  the  example). 
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A  synthetic-aperture  collection  approach  may  also  be  taken.  The  approach 
described  above  is  to  gather  N  snapshots  of  the  laser  speckle  pattern  over  a 
two-dimensional  (possibly  sparse)  array  of  detectors.  Using  a  synthetic-aperture 
approach,  only  a  one-dimensional  array  of  detectors  is  required.  For  simplicity, 
suppose  that  the  1-D  array  is  along  the  ya  axis  and  that  the  speckles  move  in 
the  xa  direction  with  velocity  va  =  2 u>R  ( =  2000  km/s  for  the  example).  At 
any  instant  in  time  t,  the  array  is  collecting  a  1-D  cut,  |Fn(oaf,ya)|  >  t«rough 
the  speckle  pattern.  As  the  speckle  pattern  sweeps  across  the  1-D  array,  ad¬ 
ditional  slices  of  speckle  pattern  fall  on  the  detector  array.  For  the  case  of  cw 
illumination,  the  detectors  must  read  out  each  sample  of  the  intensity  in  a 
time  Ah  <  X/(8i oDa)  (  =  25  ns  for  the  example).  Then  the  total  data  rate  is 
( DJdaVMi  (=  2  x  109  samples  per  second  for  the  example).  After  time  Ma 
=  DJva  (=  2.5  (is  for  the  example),  a  two-dimensional  array  of  values  of 
\Fn(xa,ya)\2  has  been  read  out,  equivalent  to  a  square  area  of  length  Da  After 
time  NMa,  a  rectangular  array  of  values  of  | F n(Xo,y0)|  has  been  read  out  o 
size  NDa  *  Da-  This  rectangular  array  can  be  subdivided  into  N  square  arrays 
of  size  Da  x  Da,  each  equivalent  to  one  of  the  snapshots  (realizations)  of  the 
speckle  pattern  that  would  have  been  measured  if  we  had  a  filled  2-D  array. 
Even  better  than  dividing  the  rectangle  into  squares  would  be  to  compute  the 
covariance  function  C(u,v )  by  summing  over  the  products  of  all  pairs  of  inten¬ 
sity  measurements  separated  by  the  distance  (u,v).  This  would  increase  the 
SNR  somewhat  for  spatial  frequencies  having  components  in  the  u  direction, 
since  separations  would  be  used  that  would  otherwise  have  crossed  the  bound¬ 
aries  between  the  squares  of  data  forming  the  snapshots.  This  idea  of  synthe¬ 
sizing  a  2-D  array  from  a  1-D  array  has  the  advantage  of  requiring  many 
fewer  detectors  than  a  real  2-D  array,  but  the  detectors  must  be  read  out  very 
fast  and  must  be  oriented  appropriately  with  respect  to  the  direction  of  the 

speckle  velocity. 


1.5.4. 2  Transmitter  Requirements.  The  received  illumination^  assumed 
to  be  quasimonochromatic,  i.e.,  the  intensity  pattern  |F„(x0,y0)|  should  be 
constant  over  the  wavelength  band  in  use.  The  requirements  on  the  coherence 
length  of  the  laser  almost  always  supersede  this  requirement,  which  can  thereby 
be  ignored.  The  coherence  length  required  of  the  laser  illuminator  is 


( =  10  m  for  the  example),  where  Az  is  the  depth  of  the  object  along  the  LOS. 

If  the  laser  is  pulsed,  then  the  pulse  repetition  frequency  is  required  to  be 
fast  enough  in  order  to  obtain  a  large  number  of  speckle  realizations  before 

the  object  rotates  too  much.  .  , 

The  laser  power  requirements  are  driven  by  several  factors,  including  tne 
fraction  of  the  transmitted  energy  that  falls  on  a  detector  and  the  SNR  of 
\Fi(u,u)\2  for  correlography. 


1 .5.4.3  Signol-to-Noise  Ratio.  The  SNR  of  \Fi(u,uf  (the  signal  divided  by 
the  standard  deviation  of  the^noise),  as  limited  by  the  number  of  realizations 
and  photon  noise,  is  given  by 
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SNR(u,i>)  = 
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where 

rid  =  the  average  number  of  photons  per  detector 

M  =  the  number  of  detectors  (pixels)  per  aperture-plane  speckle  area 

Mrid  -  the  average  number  of  photons  per  speckle 

ns  =  the  number  of  speckles  within  the  area  of  the  aperture. 

In  this  formula,  J|x|2  is  a  function  of  ( u,v ).  Figure  1.44  shows  two  families  of 
plots  of  SNR(i/,iO  as  a  function  of  light  level  for  two  values  of  |pt|2  (0.25  and 
O.OD4  and  several  values  of  Mrid,  the  number  of  photons  per  speckle.  The  first 
three  terms  in  the  denominator  of  the  expression  for  SNR(«,u)  are  due  to  the 
approximation  of  the  ensemble  average  for  the  autocovariance  by  a  finite  sum 
over  N  realizations,  and  the  last  two  terms  are  due  to  photon  noise.  For  |p.|2 
«  1,  and  for  M  =  4  detectors  per  speckle,  these  two  contributions  to  the  noise 
are  equal  when  rid  =  0.5  photon  per  detector  or  Mrid  =  2  photons  per  speckle. 
Consequently,  photon  noise  can  be  ignored  when  there  are  much  greater  than 
2  photons  per  speckle.  For  the  case  of  a  photon-counting  detector,  any  laser 
power  that  yields  more  than  about  8  photons  per  speckle  is  wasting  energy, 
since  at  that  point  it  is  the  number  of  realizations,  not  the  photon  noise,  that 
is  limiting  the  SNR,  as  can  be  seen  from  Fig.  1.44. 

For  the  higher  spatial  frequencies,  for  which  ||x|2  «  1,  and  for  moderate 
light  levels  (greater  than  two  photons  per  speckle),  the  SNR  expression  reduces  to 

SNR(u,o)  -  [(N/3)nsS(u,o)]vV(^)|2  .  (1.171) 

For  a  given  object  and  resolution,  this  expression  shows  the  trade-off  that  can 
be  made  between  the  number  N  of  realizations  and  the  number  nsS(u,v)  of 
redundant  measurements  of  a  given  spatial  frequency  (u,v). 


Fig.  1.44  The  SNR  versus  number  of  photons  per  speckle.  Solid  lines:  |p[2  =  0.25;  dashed 
lines:  |p.|2  =  0.01.  Number  of  realizations  (top  to  bottom):  (a)  10,000;  (b)  1000;  (c)  100;  (d)  10; 
and  (e)  1. 
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Fig  1  45  Three  types  of  objects:  (a)  four  points-a  2  x  4-pixel  rectangle  plus  three  single- 
pixel  points;  (b)  satellite — an  object  with  glints;  and  (c)  bus/RV-an  extended  diffuse  object. 


The  SNR  is  approximately  proportional  to  |m-(u,o  )|  >  which  heavily  depends 
both  on  (u,v)  and  on  the  structure  of  the  object  being  imaged.  At  one  extreme, 
an  object  consisting  of  a  single  unresolved  glint  (a  delta  function)  has  Wwl 
=  1.0  for  all  ( u,v ).  At  the  other  extreme,  a  diffuse,  constant- value  object  that 
fills  the  field  of  view  will  have  |(x(u,o)|  that  is  a  delta  function  and  is  negligibly 
small  for  (uv)  #  (0,0).  To  give  you  a  feel  for  the  values  that  |p-(w,o)|  takes  on 
for  different  types  of  objects,  Fig.  1.45  shows  three  types  of  objects  and  Fig. 
1.46  shows  |pL(ur)|,  the  spin-averaged  (azimuthally  averaged)  values  of  the 
visibility  |ji(u,u)|,  where  ur  =  («2  +  ^  )  •  A  vei7  rough  rule  of  thumb  is  that 
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where  nbp  is  the  number  of  resolved  bright  image  points  at  resolution  ur.  For 
an  extended  object  npb  is  the  number  of  2-D  resolution  elements  within  the 
area  of  the  object.  So  if,  as  the  resolution  improves  with  increasing  ur,  more 
points  in  the  diffraction-limited  image  of  an  object  are  resolved,  then  |p,(i/r)[2 
decreases  correspondingly. 

Ideally  the  laser  should  have  sufficient  power  so  that  the  number  of  detected 
photons  per  speckle  exceeds  two.  For  a  simplified  object  consisting  of  a  flat 
plate  that  is  a  Lambertian  reflector,  the  average  number  of  photons  per  speckle 
is  given  by 

2  2  2 
Mrid  =  ELi&tmt\LA{cosfti  cos0o/ir)roTpoiTopt(Tidc(a/J?  )t\q(\/hc)  ,  (1.173) 

where 

El  =  energy  of  the  laser  per  pulse  or  per  detector  integration  time 

Tatm  =  one-way  atmospheric  transmittance 

t\la  =  fraction  of  laser  intensity  falling  on  the  area  of  the  object 

0,  =  angle  between  object  surface  normal  and  illumination 

0O  =  angle  between  object  surface  normal  and  LOS 

rQ  =  average  reflectivity  of  object 

Tpoi  =  fraction  of  light  reflected  in  desired  polarization 

Topt  =  transmittance  of  receiver  optics  (if  any) 

r\dd2a  =  area  of  one  detector 

j\q  ~  quantum  efficiency  of  the  detector 

hc/\  =  energy  of  one  photon  (Planck’s  constant  times  frequency). 

For  example,  for  the  parameters  of  Table  1.3  and  for  El  -  1  J/pulse, 
TatmTpol^opt  0.2,  T| j^A  —  0.1,  0(  =  0o  =  0,  To  =  0.5,  T\dda  —  0.005  m  ,  T|q  — 
0.1,  and  hc/K  =  2  x  10' 19  J/photon,  we  have  Mns  =  8  photons  per  speckle 
per  pulse,  which  is  in  the  high  light-level  regime. 

If  the  receiver  is  ground  based,  then  atmospheric  scintillation  could  be  a 
limiting  factor,  since  it  directly  changes  the  measured  aperture-plane  inten¬ 
sity.  However,  since  many  speckle  patterns  are  averaged,  the  effects  of  scin¬ 
tillation  tend  to  cancel  and  are  not  expected  to  be  a  limiting  factor  for  imaging 
objects  directly  overhead. 

1.5.5  Conclusion 

The  use  of  laser  illumination  opens  the  possibility  for  many  novel  ways  to 
image  objects  with  resolution  one  to  two  orders  of  magnitude  finer  than  con¬ 
ventional  imaging  systems  that  are  limited  by  the  blurring  effects  due  to 
atmospheric  turbulence.  These  imaging  methods,  which  spring  from  the  her¬ 
itage  of  holography,  are  especially  useful  for  obtaining  fine-resolution  images 
of  earth-orbiting  satellites.  These  same  methods  can  also  be  applied  to  other 
imaging  scenarios,  ranging  from  microscopy  to  imaging  of  the  ground  from 
space.  While  these  methods  have  been  investigated  and  demonstrated  to  vary¬ 
ing  degrees  of  sophistication,  they  are  just  now  on  the  verge  of  becoming 
operationally  feasible. 
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ABSTRACT 

We  describe  two  new  methods  for  phasing  arrays  of  heterodyne  receivers.  Both  can  be  used  when  the 
arrays  are  sparse  and  distributed.  One  is  based  on  the  iterative  transform  algorithm  using  a  support 
constraint  and  the  other  on  maximizing  image  sharpness.  Both  work  well:  they  require  a  modest  number  of 
speckle  realizations  with  the  same  aberrations  and  are  relatively  immune  to  measurement  noise. 

1  INTRODUCTION 


If  we  measure  an  optical  field,  scattered  from  a  laser-illuminated  object,  with  an  array  of  heterodyne 

|  • _ .r  j _ Ut?  norfAtmin nr  a  cintPmft  TlTOTIfl 
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receivers/  then  we  can  reconstruct  an  image  of  the  object  by  performing  a  suitable  propagation 
transformation  (typically  a  Fresnel  or  Fourier  transform)  in  the  computer.2  However,  if  phase  errors  ;  are 
present  in  the  measurements,  then  a  blurred  image  will  result.  Phase  errors  can  be  induced  by  many th mgs 
including  atmospheric  turbulence  and  path-length  errors  between  the  field  from  the  object  and  the  local 
oscillator  within  the  receiver.  Approaches  to  phasing  the  array  without  a  beacon  either  require  the  array  of 
heterodyne  receivers  to  be  on  a  filled  regular  grid,3’4  or  to  have  a  special  pattern  or  alternatively  require  a 
low-resolution  image  be  available.5 

In  this  paper  we  describe  two  new  algorithms  for  phasing  heterodyne  arrays  that  work  well  for  sparse 
distributed  arrays  (and  will  work  for  filled  arrays  as  well),  requires  no  low-resolution  image,  and  is 
relatively  insensitive  to  noise.  One  is  based  on  the  iterative  transform  algorithm  using  a  support 
constraint6’7  and  the  other  on  maximizing  image  sharpness.8- 

9.  IMAGING  MODEL 


We  measure  the  optical  fields  (ignoring  measurement  noise) 

Gdk(u)  =  Fk(u)  exp[i<j>e(u)] 


(1) 


for  k  =  1  ,  K,  where  Fw  is  the  ideal  complex  field  (without  phase  errors)  for  realization  k,  4>e  is  the 

phase  error,  and  u  is  a  2-D  coordinate  (the  pixel  index)  in  the  measurement  plane.  Note  that  we  assume  the 
phase  error  is  the  same  for  all  K  realizations.  Different  speckle  realizations  are  obtained  if  the  object  rotates 
slightly  or  translates.  If  we  estimate  the  phase  error  as  <)(u),  then  our  estimates  of  the  optical  fields  are 


Gv(u)  =  Gdk(u)  exp[-i<{>(u)] 


(2) 


our  coherent  images  are 
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(4) 


gk(x)  =  jF_1[Gk(u)]  =  N-2£Gk(u)  exp(i27iu-x/N)  ,  (3) 

u 

where  irx  is  a  vector  dot  (inner)  product,  and  our  speckle-reduced  averaged  image  is 

gj(x)  =  K^XlgkOOl2  * 
k 

where  summation  over  k  is  for  k  =  1, K. 

3.  JTF.RATIVF.  TRANSFORM  ALGORITHM 

Our  first  method  to  estimate  the  phase  error  is  a  generalization  of  the  iterative  transform  algorithm.6  In 
summary,  one  iteration  of  the  generalized  algorithm  consists  of  the  following  steps: 

(1)  Calculate  the  Fourier  estimates:  Gk(u)  =  G^u)  exp[-i(j)(u)]  (for  k  =  1  ...  K). 

(2)  Inverse  FFT  them  to  compute  the  gk(x). 

(3)  Apply  the  support  constraint:  gsk(x)  =  s(x)  gk(x). 

(4)  FFT  to  compute  the  Gsk(u). 

(5)  Compute  a  new  phase-error  estimate:  <{>(u)  =  arg{^  Gdk(u)  Gsk(u)  J. 

k 

For  the  initial  iteration  the  phase  error  estimate  is  zero.  The  support  constraint  can  be  known  a  priori 
estimated  from  the  autocorrelation  support  (which  can  be  calculated  from  the  given  data),  *  or  estimated 
from  gT(x).  To  estimate  the  support  from  gl(x),  the  current  estimate  of  the  averaged  image  we  could 
smooth  it,  giving  s(x)  =  w(x)  *  gl(x),  and  possibly  threshold  it.  The  smoothing  function  w(x) cou id I  e 
for  example,  a  uniform  disk  having  a  radius  of  a  few  pixels.  The  new  phase-error  estimate,  4>(u),  given  m 
step  5  is  like  a  weighted  average  of  the  difference  in  phase  between  Gdk(u),  the  measured  data,  and  Gsk(u), 
the  current  estimate  of  optical  fields  constrained  by  the  image  support  constraint. 

Steps  (2)  to  (4)  can  also  be  performed  by  a  direct  convolution  in  the  measurement  plane. 


Gsk(u)  =  Gk(u)  *  S(u)  =  £  S(u2)  Gk(u  -  uj)  , 


(5) 


where  S(u)  =  7[ s(x)],  which  would  be  faster  to  compute  for  the  case  of  sparse  arrays. 

d  maximizing  sharpness 


The  second  approach  is  to  find  the  phase-error  estimate  <J)(u)  which  maximize  the  sharpness, 

S  1  =  Z  [g^x)]2  =  LtXlgkWI2]2  ’  (6) 

X  X  k 

of  the  image  estimate.  Alternatively,  we  can  maximize  S1W,  a  weighted  version  of  the  sharpness,  based  on 
a  filtered  version  of  the  averaged  image, 

giw(x)  =  y-1(wW7feiW]l  (7) 
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where  W(u)  is  a  weighting  function  (for  example,  a  Wiener  filter). 

We  maximize  the  sharpness  using  a  standard  nonlinear  optimization  technique  which  we  make 
computationally  efficient  by  using  an  analytic  expression  for  the  gradient  of  the  sharpness  with  respect  to 
the  phase-error  estimate: 


=  4N-2  Im{  V  Gk(u)  (Gk(u)  *  [GjW2])*}  ,  (8) 

3<KU)  k 


where  Gj(u)  =  5[gj(x)]. 

We  can  show  that  the  error-reduction6  version  of  our  generalized  iterative  transform  algorithm  is  similar 
in  effect  to  a  steepest-descent  gradient  search  to  maximize  sharpness  with  a  fixed  step  size.  However,  we  in 
practice  use  a  more  powerful  gradient  search  algorithm,  such  as  conjugate-gradient,  and  we  can  use  a  more 
powerful  version  of  the  iterative  transform  algorithm  as  well,  such  as  hybrid  input-output. 

5.  TMAGE  RECONSTRUCTION  EXAMPLE 

Figure  1  shows  a  computer  simulation  example  of  phasing  an  array  of  sparse,  distributed  heterodyne 
receivers.  Figure  1(a)  shows  the  object,  consisting  of  a  block  of  2x3  pixels  and  three  separated  points, 
embedded  in  a  64x64  array.  We  created  complex- valued  speckled  images  by  replacing  each  pixel  in  the 
object  by  a  circular  complex  Gaussian  random  number  with  variance  equal  to  the  square  root  of  the  original 
pixel  value;  Fourier  transforming,  multiplying  by  the  a  sparse  aperture,  and  inverse  transforming  to 
produce  the  image.  The  sparse  aperture  we  used  for  this  experiment  has  21  pixels  within  a  8x8  area,  which 
is  shown  in  Figure  1(b),  magnified  by  a  factor  of  four  (only  the  central  16x16  pixels  in  its  64x64  array  are 
shown).  This  sub-aperture  pattern  was  designed  by  a  genetic-based  algorithm  by  P.  Henshaw.  An 
average  of  ten  such  diffraction-limited  speckled  images  is  shown  in  Figure  1(c).  Figure  1(d)  shows  a 
Wiener-Helstrom1 3  filtered  version  of  the  averaged  image,  from  which  one  can  discern  the  major  parts  of 
the  image.  Figure  1(e)  shows  the  blurred  averaged  image  that  resulted  when  we  added,  to  the  Fourier 
transform,  random  phase  errors  that  were  uniformly  distributed  Gaussian  random  numbers  with  standard 
deviation  %  radians.  Figure  1(f)  shows  the  averaged  image  corrected  by  the  image  sharpening  algorithm 
and  Figure  1(g)  shows  a  Wiener-filtered  version  of  the  corrected  image.  The  corrected  image  is  clearly 
superior  to  the  blurred  image  and  strongly  resembles  the  diffraction-limited  image.  We  obtained  comparable 
results  using  the  generalized  iterative  transform  algorithm.  The  algorithms  usually  produced  images  with 
better  quality  as  the  number  of  realizations  increases.  It  was  unreliable  when  using  only  one  to  four 
realizations.  Tests  with  noise  showed  that  the  algorithm  is  not  very  sensitive  to  noise. 

6  CONCLUSIONS 

We  have  derived  two  new  algorithms  for  phasing  an  array  of  heterodyne  receivers,  and  have 
demonstrated  that  they  can  work  well  for  sparse  distributed  arrays,  require  a  modest  number  of  speckle 
realizations,  and  are  robust  in  the  presence  of  measurement  noise. 
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(a) 


Figure  1.  Computer  Simulation  of  Sparse  Heterodyne  Array  Phasing. 

(a)  Object;  (b)  21 -element  sparse,  distributed  array  (within  an  8  x  8  -sample  area); 
(c)  average  of  ten  diffraction-limited  speckled  images;  (d)  Wiener-filtered  version 
of  (c);  (e)  average  of  ten  aberrated  speckled  images;  (f)  average  of  ten  speckled 
images  corrected  by  image-sharpening  phase-up  algorithm;  (g)  Wiener-filtered 
version  of  (f). 
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ABSTRACT 

Three-dimensional  imaging  provides  profile  information  not  available  with  conventional  two- 
dimensional  imaging.  Many  three-dimensional  objects  of  interest  are  opaque  to  the  illuminating 
radiation,  meaning  that  the  object  exhibits  surface,  as  opposed  to  volume,  scattering.  We  inves¬ 
tigate  the  use  of  an  opacity  constraint  to  perform  three-dimensional  phase  retrieval.  The  use  of 
an  opacity  constraint  in  conjunction  with  frequency-diverse  pupil-plane  speckle  measurements  to 
reconstruct  a  three-dimensional  object  constitutes  a  novel  unconventional-imaging  concept.  This 
imaging  modality  avoids  the  difficulties  associated  with  making  phase  measurements  at  a'  cost  of 
increased  computations. 


1.  INTRODUCTION 

Marron  and  Schroeder  [1,2]  have  demonstrated  that  three-dimensional  imaging  can  be  accom¬ 
plished  by  sequentially  illuminating  an  object  with  different  laser  frequencies  and  measuring  the 
far-field  speckle  pattern  for  each  of  the  illuminating  frequencies  [1].  In-phase  and  quadrature  mea¬ 
surements  provide  optical-field  values  for  the  laser-speckle  patterns.  Properly  formatted,  these  data 
form  a  three-dimensional  Fourier-volume  (or  Fourier-aperture)  representation  of  the  illuminated 
object.  A  simple  three-dimensional  Discrete  Fourier  Transform  (DFT)  can  be  performed  to  provide 
a  three-dimensional  representation  of  the  object.  This  three-dimensional  lensless-imaging  concept 
is  referred  to  as  Holographic  Laser  Radar  (HLR). 

The  collection  of  HLR  data  could  be  considerably  simplified  if  instead  of  collecting  field  mea¬ 
surements  for  the  speckle  images,  intensity  measurements  were  collected.  Intensity  measurements 
are  straightforward  and  eliminate  the  need  for  interference  with  a  reference  beam  with  precision 
alignment,  tracking,  and  phase  stability.  A  multiple-frequency  speckle-intensity  data  set  would  pro¬ 
vide  three-dimensional  Fourier-magnitude  information  or,  equivalently,  the  three-dimensional  object 
autocorrelation.  In  order  to  recover  a  literal  three-dimensional  image,  a  phase-retrieval  algorithm 
is  required.  Therefore,  this  proposed  imaging  concept  trades  complexity  and  cost  in  hardware  for 
increased  computing. 

Phase  retrieval  requires  some  type  of  a  priori  information  about  the  object.  Two-dimensional 
complex- valued  objects  have  been  recovered  using  phase  retrieval  with  a  support  constraint  [3], 
although  this  is  a  challenging  problem.  It  is  well  known  that  the  uniqueness  properties  of  two- 
dimensional  phase  retrieval  are  much  better  than  for  the  one-dimensional  problem  [4].  We  conjecture 
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that  three-dimensional  phase  retrieval  with  a  support  constraint  is  better  conditioned  than  its  two- 
dimensional  counterpart. 

An  additional  constraint  that  has  great  promise  in  the  three-dimensional  imaging  case  is  an 
opacity  constraint.  An  opaque  object  is  one  that  exhibits  only  surface  scattering  and  no  volume 
scattering  (over  volumes  that  extend  beyond  the  desired  range  resolution).  The  reflectivity  function 
for  an  opaque  object  is  confined  to  a  two-dimensional  curved  (possibly  discontinuous)  surface, 
embedded  in  a  three-dimensional  space.  In  this  paper  we  explore  the  use  of  an  opacity  constraint  to 
perform  phase  retrieval.  We  have  previously  used  an  opacity  constraint  to  perform  superresolution 
[5]. 


The  opacity  constraint  is  a  special  type  of  support  constraint.  It  is  a  “quality  of  support” 
constraint  -  the  actual  location  of  the  support  is  not  given,  although  the  object  is  known  to  be 
confined  to  a  two-dimensional  curved  surface.  This  constraint  promises  to  be  very  powerful  since 
it  greatly  reduces  the  class  of  feasible  objects  from  which  to  choose  an  estimate.  Moreover,  there 
are  many  imaging  applications  in  which  the  objects  will  be  known  with  confidence  to  be  opaque. 
For  example,  space  objects,  balistic  missiles,  aircraft,  and  a  multitude  of  industrial-inspection  parts 
qualify  as  opaque  objects.  Most  objects  in  our  everyday  experience  satisfy  the  opaque  condition. 
The  constraint  is  invalid  for  objects  with  distributed  volume  scatterers  such  as  translucent  or  fog- like 
objects. 

The  use  of  an  opacity  constraint  in  conjunction  with  frequency-diverse  pupil-plane  speckle  mea¬ 
surements  to  reconstruct  a  three-dimensional  object  constitutes  a  novel  unconventional-imaging 
concept. 


2.  REPRESENTATION  OF  OPAQUE  OBJECTS 


Consider  a  three-dimensional  opaque  object,  f(x,  y ,  z),  defined  on  an  object-centered  coordinate 
system.  Let  x  and  y  be  the  cross-range  (angle-angle)  coordinates  and  let  z  be  the  range  coordinate 
that  is  co-aligned  with  the  illumination  direction.  Since  the  collected  data  will  be  sampled  and 
since  object  reconstructions  must  be  performed  with  a  digital  computer,  we  will  adopt  a  discrete 
representation  for  the  object.  Accordingly,  we  require  the  coordinates  ( x,y,z )  to  take  on  integer 
values.  Because  of  the  opaque  nature  of  the  object,  only  radiation  from  reflecting  sources  in  the 
illuminated  surface  will  contribute  to  the  received  signal.  Hidden  surfaces  do  not  contribute.  Let 
h(x,  y)  denote  the  height  of  the  object  in  the  ^-dimension.  Because  the  object  is  confined  to  a 
two-dimensional  curved  surface  embedded  in  three-dimensional  space,  it  can  be  represented  with 
delta-function  notation  [5,6], 

f{x,y,z)  =  r(x,y)5[z- h(x,y)\  ,  x,  y,  z  €  {0, 1, 2, . . .  N  -  1}  ,  (1) 


where  r(x,  y)  is  the  complex  surface  reflectivity,  N  is  the  total  number  of  samples  in  each  dimension, 
and  <$(•)  is  used  to  represent  the  Kronecker  delta, 


S(z  -  z0 )  = 


1,  z  =  z0 
0,  z^  z0 


(2) 
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Recall  that  the  collected  Fourier  intensity  data  can  be  transformed  to  compute  a  three-dimensional 
object  autocorrelation.  In  order  to  avoid  aliasing,  we  require  the  object  to  have  finite  support  with 
sufficient  zero  padding, 

f{xty,z)  =  0,  x,y,z€  y  +  l,...,N-lj  .  (3) 

The  rectilinear  object-support  bound  given  in  Eq.  (3)  is  not  likely  to  be  a  “tight”  bound.  More 
restrictive  or  tighter  bounds  on  object  support  can  be  determined  from  the  known  object  auto¬ 
correlation  by  using  a  generalization  of  methods  developed  for  two-dimensional  phase  retrieval  [7]. 
Any  three-dimensional  support  bound  can  be  expressed  in  terms  of  a  two-dimensional  angle-angle 
support  bound  and  a  one-dimensional  range  support  bound  that  depends  upon  angle-angle  position. 

3.  STATEMENT  OF  PROBLEM 

The  Fourier  representation  of  the  object  is  found  with  the  DFT, 

N- 1  N-1N-1 

F(u,v,w)  —  11,  5E  f(x,y,z)exp{-i2n(ux  +  vy  +  wz)/N}  (4) 

a;=0  y= 0  z=0 

N—l  N—l  N—l 

=  EEE^:#  —  h(x,  y )]  exp{— i2Tr(ux  +  vy  +  wz)/N}  (5) 

x=0  y—0  z= 0 

N- 1  N-l 

=  £2>U  ,  y)  exp{— i2-K[ux  +  vy  +  wh(x,  y)]/N}  ,  (6) 

z=0  y=0 

where  we  have  used  the  Kronecker  delta  to  eliminate  the  summation  over  z.  Since  we  detect 
intensities  in  the  Fourier  domain,  a  noiseless  measurement  would  be  given  by  the  squared  modulus 
of  the  DFT  of  the  object,  \F(u,  v,  io)|2.  The  actual  detected  data  will  be  corrupted  by  noise.  The 
detected  data  are  represented  by 

D(u,v,w)  =  N {|F(u, u, w)\2}  ,  (7) 

where  the  noise  operator  Af{-}  corrupts  the  argument  according  to  an  appropriate  noise  model. 
Additive  Gaussian  or  Poisson  noise  models  are  appropriate  when  detector  readout  noise  or  photon 
noise,  respectively,  are  the  dominant  noise  sources. 

The  problem  that  we  wish  to  address  can  now  be  stated.  Given  the  data,  D(u,v,w),  estimate 
the  object  height  function,  h(x,y ),  and  the  object  surface  reflectivity  function,  r(x,y),  over  the 
angle-angle  object-support  bound.  Performing  such  an  estimate  implicitly  provides  an  estimate  for 
the  object’s  Fourier  phase.  Therefore  this  problem  may  be  properly  viewed  as  a  phase-retrieval 
problem. 


4.  MAXIMUM-LIKELIHOOD  ESTIMATION 

When  detector  readout  noise  dominates,  the  noise  operator  has  the  effect  of  adding  a  realization 
of  a  Gaussian  noise  process, 

D(u,v,w )  =  |F(u,u,  w)\2  +  n(u,v,w)  .  (8) 
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Assuming  a  common  noise  variance  and  statistical  independence  between  differing  samples  of  the 
noise  process,  it  is  straightforward  to  express  the  probability  density  function  for  the  recorded  data 
[8]  .  The  associated  log-likelihood  function  is  readily  found  to  be  [8] 


N-1N-1N-1 

L[r(x,y),h(x,y )]  =  -  EEE  [D(u,V,w)  -  lP{u,V,w)\\ 


N-1N-1N-1 

N-1N-1 

2" 

EEE 

D(u,  v,  w)  — 

X  X  r{x,  y)  exp{-i2n[ux  +  vy  +  wh(x,  y)]/N] 

u=0  v=0  w= 0 

O 

II 

o 

II 

w 

(9) 

(10) 


where  the  caret  symbol  indicates  an  estimated  quantity.  Recall  that  the  object  surface  reflectivity  is 
complex  valued.  It  is  convenient  to  explicitly  denote  the  real  and  imaginary  parts  of  the  reflectivity, 


r(x,y)  =  rr(x,y)  +  ir{(x,y)  .  (11) 

We  can  use  standard  nonlinear  optimization  methods  to  search  for  rr(x,y),  fi(x,y),  and  h(x,y ) 
((x,y)  £  the  angle-angle  object-support  bound)  that  maximize  the  objective  function,  L ,  to  yield  a 
maximum-likelihood  estimate.  If  the  angle-angle  object  support  bound  is  comprised  of  N/2  by  N/2 
samples,  then  L  is  defined  on  a  parameter  space  of  dimension  3fV2/4.  The  range  support  bound 
that  derives  from  the  three-dimensional  support  bound  could  be  used  to  constrain  the  search  over 
h(x,y). 


There  are  a  variety  of  nonlinear  optimization  algorithms  that  could  be  used  to  maximize  the 
log-likelihood  function.  One  choice  is  the  well-known  conjugate-gradients  algorithm  [9].  As  its 
name  suggests,  this  algorithm  makes  use  of  the  gradient  of  the  objective  function.  We  have  derived 
closed-form  expressions  for  the  partial  derivatives  that  constitute  the  gradient.  These  expressions 
afford  a  gradient  computation  that  is  more  accurate  and  computationally  much  more  efficient  than 
would  be  rendered  using  the  method  of  finite  differences.  The  partial  derivatives  of  the  log-likelihood 
function  with  respect  to  the  real  and  imaginary  parts  of  the  surface  reflectivity  and  with  respect  to 
the  height  function,  at  the  angle-angle  location  ( x0,ya ),  are  given  by 

dL 

drr{x0,y0) 

( N-1N-1N-1  ~  ^  'i 

4Re<  X  X  X  [d{u,v,w)~  |F(u,u,u;)|2]  F*(u,v,w)  exp{-i2n[ux0  +  vy0  +  wh(x0,y0)]/N}  >,  (12) 
l  u= 0  v=0  w=0  J 

dL 

dfi{x0,  y0) 

( N-1N-1N-1  _  'I 

— 41m  <  X]  \D{u,v,w)  -  \F(u,v,w)\2^  F*(u,v,w)  exp{-i2Tr[ux0  +  vy0  +  wh(x0,y0)]/N}> ,  (13) 

l  u~ 0  v—0  w=0  J 
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dL 


8i T 


dh(xo,y0) 


N-lN-lN-l 

EEE 

u=0  v=0  w—0 

F*(u ,  v,  w)  wexp{-i2TT[ux0  +  vy0  +  wh(x0,  J/0)]/IV}l  ,  (14) 


^Im  lr(xo,y0)Y^  Y  Y  D(u,v,w)  -  \F(u,v,w)\2. 


where  the  operators  Re{-}  and  Im{-}  take  the  real  and  imaginary  parts  of  the  argument,  respectively. 
Notice  that  these  partial-derivative  expressions  each  have  the  form  of  a  two-dimensional  DFT  (over 
u  and  v)  followed  by  a  summation  over  w  for  each  (x0,  y0).  This  affords  efficient  evaluation  the  entire 
gradient  using  the  Fast  Fourier  Transform  (FFT).  The  three-dimensional  function  F(u,v,w)  must 
be  computed  to  evaluate  the  objective  function  and  the  gradient.  We  find  that  N  two-dimensional 
FFTs  are  needed  to  compute  F(u,v,w).  Careful  examination  of  Eqs.  (12)  — (14)  reveals  that  the 
entire  gradient  can  be  computed  with  an  additional  N  two-dimensional  FFTs  followed  by  2SN 2 
summations  over  w ,  where  S  is  the  proportion  of  pixels  in  the  angle-angle  object  support  relative 
to  the  number  of  samples  in  the  total  angular  field  of  view.  In  the  case  of  the  rectilinear  support 
expressed  in  Eq.  (3),  S  =  1/4.  Of  course  additional  computational  overhead  is  required  to  perform 
complex  multiplies  and  sums.  Closed-form  expressions  for  the  log-likelihood  and  its  gradient  can 
also  be  derived  for  the  case  of  photon  noise. 

5.  SIMULATIONS 

We  have  performed  two  simulations  that  demonstrate  three-dimensional  phase  retrieval  using 
an  opacity  constraint.  For  both  of  the  simulations  we  chose  N  =  32  and  a  maximum  object  extent 
of  16  samples  in  each  dimension,  as  required  by  Eq.  (3).  This  guarantees  that  the  Fourier-intensity 
data  are  adequately  sampled.  In  addition,  for  this  initial  demonstration  we  used  noiseless  Fourier- 
intensity  data,  giving  n(u,v,w)  =  0  in  Eq.  (8). 

The  first  simulation  utilized  a  simple  object  consisting  of  6  separated  points  of  equal  reflect¬ 
ing  strength.  The  angle-angle  (x-y)  view  of  the  magnitude  of  the  surface  reflectivity  is  shown  in 
Figure  1(A).  The  phases  of  these  points  are  independent,  identically-distributed  samples  from  a 
probability  distribution  function  (PDF)  that  is  uniformly  distributed  on  the  interval  [0,27 r].  The 
height  distribution  and  reflectivity  of  these  points  are  summarized  in  Table  1.  Coordinates  are 
referenced  from  (0, 0)  in  the  upper  left  corner.  Despite  the  discrete  representation  of  the  object  in 
the  x  and  y  dimensions,  the  ^-dimension  heights,  h(x,y),  are  not  restricted  to  integer  values. 

The  phase-retrieval  problem  consists  of  estimating  a  set  of  real-valued  quantities,  rr(x,y), 
ri(x,y ),  and  h(x,y),  defined  on  a  discrete  ( x,y )  coordinate  system.  The  conjugate-gradients  al¬ 
gorithm,  used  to  estimate  the  object  parameters,  is  initiated  with  a  first  guess  of  the  unknown 
parameters.  We  began  with  an  initial  estimate  of  h(x,y )  =  0  and  independent,  identically  dis¬ 
tributed  samples  from  a  PDF  that  is  uniform  on  the  interval  [0,1]  for  rr(x,y)  and  ri(x,y).  The 
magnitude  and  phase  of  the  initial  reflectivity  guess  are  shown  in  Figures  1(B)  and  1(C),  respec¬ 
tively.  The  magnitude  of  the  initial  reflectivity  guess  ranges  between  0.083  and  1.34,  and  the  phase 
varies  from  0.001  to  1.567  radians.  By  restricting  the  initial  guess  of  rr(x,  y)  and  r,-(x,  y )  to  positive 
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(C) 


(D) 


Figure  1.  6-point  object  simulation.  (A)  Angle-angle  ( x-y )  view  of  the  magnitude  of  the 
reflectivity  of  a  6-point,  complex-valued  object;  (B)  Magnitude  of  the  initial 
reflectivity  estimate;  (C)  Phase  of  the  initial  reflectivity  estimate;  (D)  Recon¬ 
structed  object. 
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(x,y) 

rr{x,y) 

rt{x,y) 

h(x,y) 

(1,3) 

-0.976 

0.216 

0.0 

(2,7) 

0.411 

0.911 

0.3 

(7,8) 

0.994 

0.111 

0.9 

(3,12) 

-0.988 

0.153 

1.5 

(8,14) 

0.745 

-0.666 

2.1 

(11,12) 

0.089 

0.996 

3.0 

Table  1:  6-Point  Object  Description. 


values,  we  have  reduced  the  chance  of  biasing  the  simulation  by  starting  too  close  to  the  true  solu¬ 
tion.  No  a  priori  knowledge  of  the  support  of  the  unknown  object  was  used  during  the  estimation 
of  its  reflectivity  and  height  distribution,  so  a  total  of  16  x  16  x  3  =  768  parameters  were  estimated 
simultaneously. 

The  angle-angle  view  of  the  estimated  reflectivity  magnitude  after  about  400  iterations  of  the 
conjugate-gradients  algorithm  is  shown  in  Figure  1(D).  The  reconstruction  is  virtually  perfect,  and 
is  limited  only  by  computational  precision.  There  are  several  features  of  the  reconstruction  to  note. 
First,  the  reconstructed  object  is  shifted  from  its  true  location.  A  feature  of  all  phase  retrieval 
algorithms  is  a  translation  ambiguity  among  the  set  of  valid  object  solutions.  A  translation  of  the 
object  adds  a  linear  phase  term  to  the  complex- valued  Fourier  transform,  but  does  not  affect  the 
Fourier  intensity.  There  is  also  a  180-degree  object-rotation  ambiguity.  Such  a  rotation  gives  rise 
to  a  Fourier  transform  that  is  the  complex-conjugate  of  the  original  Fourier  transform,  but  which 
leaves  the  Fourier  intensity  unchanged. 

The  height  estimates  were  also  virtually  perfect.  It  is  of  interest  that  the  non-integer  height 
values  were  estimated  very  accurately,  indicating  that  in  addition  to  phase  retrieval,  we  were  ac¬ 
complishing  superresolution  in  range.  Not  only  is  the  object  estimate  translated  in  ( x,y ),  but  it  is 
also  translated  in  the  z  direction  by  a  constant  value  of  —1.2.  Also,  whereas  the  magnitude  of  the 
estimated  reflectivity  matches  the  true  reflectivity,  the  phase  of  the  estimate  differs  by  the  constant 
value  of  0.64  radians  at  each  point.  This  is  yet  another  ambiguity  associated  with  the  phase  retrieval 
of  complex-valued  objects.  A  constant  phase  bias  across  the  object,  however,  is  typically  irrelevant. 

For  the  case  of  noiseless  Fourier-intensity  measurements  we  have  shown  that  an  array  of  points 
with  a  modest  height  distribution  can  be  reconstructed  perfectly  with  no  a  priori  knowledge  about 
the  support  or  height  distribution  of  the  points,  except  for  the  use  of  an  opacity  constraint.  We  have, 
however,  found  that  the  conjugate-gradients  algorithm  is  subject  to  stagnation  in  what  we  believe 
to  be  local  extrema  in  the  log-likelihood  function,  and  a  successful  reconstruction  was  obtained  only 
after  several  random  initial  guesses  for  the  reflectivity.  Also,  we  have  observed  that  the  relative 
distribution  of  the  points  in  the  z-dimension  has  a  bearing  on  the  success  of  the  reconstruction.  We 
chose  a  distribution  over  a  range  of  3  samples,  and  we  have  noted  that  as  this  range  is  expanded, 
stagnation  becomes  more  likely  when  using  a  naive  initial  guess  of  h(x,y)  =  0.  Thus,  we  believe 
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that  the  development  of  methods  that  yield  improved  initial  estimates  and  support  bounds  that 
will  constrain  the  range  of  height  parameters  will  be  important  for  successful  reconstructions  for 
surfaces  that  have  large  variation  in  depth.  Our  second  simulation  provides  an  example  of  how  such 
a  height  constraint  could  be  derived  and  utilized. 

The  second  simulation  example  involves  an  object  that  is  a  curved  surface  with  an  angle-angle 
(. x-y )  support  that  is  a  5  x  5-pixel  right  triangle.  To  simplify  the  height  distribution,  we  made  the 
height  of  the  object  a  function  of  the  ^-dimension  alone.  The  support  of  the  object  projected  in 
each  of  the  3  directions  is  shown  in  Figure  2.  When  the  object  is  projected  in  the  redirection,  we 
obtain  a  single  pixel  of  ^-direction  support  for  each  position  in  the  x  dimension,  as  illustrated  in 
Figure  2(A).  We  conclude  from  this  that  the  height  is  constant  in  y  for  each  x  and  is  symmetric 
about  the  x  =  8  plane.  The  height  of  the  object  varies  over  a  5-pixel  interval,  with  a  2-pixel  height 
gap  between  the  ^-dimension  edges  and  their  neighbors. 

If  no  support  constraint  is  used  along  with  a  naive  initial  guess  of  h(x,y)  —  0,  the  parameter 
estimation  stagnates  and  the  object  cannot  be  reconstructed.  Thus,  the  reconstruction  of  this 
more  complicated  object  might  benefit  from  a  support  constraint.  Techniques  for  bounding  the 
object  support  using  the  autocorrelation  support  have  been  developed  for  two-dimensional  phase- 
retrieval  problems  [7].  We  can  apply  these  same  techniques  to  the  support  of  the  three-dimensional 
autocorrelation  projected  into  2  dimensions  in  each  of  the  3  directions.  Unfortunately,  in  doing 
so  we  lose  registration  information  between  the  3  resulting  2-D  support  constraints  and,  thus, 
cannot  merge  these  results  into  a  tight  3-D  support  constraint.  Flowever,  for  this  particular  object 
for  which  the  height  is  constant  along  the  y-direction,  the  triple-intersection  rule  [7]  applied  to 
the  support  of  the  y-projected  autocorrelation  yields  the  exact  ^-support  shown  in  Figure  2(A). 
This  somewhat  surprising  result  suggests  that  opacity  could  be  extremely  useful  in  finding  support 
bounds.  The  derived  height  support  in  turn  yields  the  true  solution  for  the  height  and  indicates 
that  the  ^-dimension  support  is  5-pixels  wide.  Merging  this  information  with  the  support  of  the 
autocorrelation  projected  along  the  ^-direction  into  the  x-y  plane  yields  the  x-y  object  support 
bound  shown  in  Figure  3(A).  This  bound  could  further  be  reduced  using  a  two-point  intersection 
rule,  but  we  found  that  this  was  not  necessary. 

The  magnitude  of  the  reflectivity  of  the  true  object  is  shown  in  Figure  3(B).  The  initial  guess 
for  the  reflectivity  is  found  in  Figure  3(C),  and  the  initial  guess  for  the  height  was  set  to  the  values 
yielded  from  the  x-z  support  derived  from  the  triple-intersection  of  the  projected  autocorrelation. 
The  initial  guess  of  the  height  was  set  constant  in  the  y-direction  over  the  support  defined  in 
Figure  3(A).  The  x-y  support  constraint  reduces  the  number  of  parameters  to  be  estimated  to 
114.  After  approximately  800  iterations  of  the  conjugate-gradients  algorithm,  the  estimates  were 
perfect  to  within  the  limits  of  computational  precision.  Despite  our  knowledge  of  the  height  from 
the  support  estimate,  we  still  allowed  the  height  estimates  to  vary  during  minimization.  The  final 
estimates  actually  moved  a  fixed  distance  from  their  initial  settings,  but  remained  constant  in  the 
y-direction.  The  final  estimate  of  the  reflectivity  is  shown  in  Figure  3(D),  and  is  essentially  identical 
to  the  true  reflectivity  in  Figure  3(B). 

Although  it  is  unlikely  that  the  height  of  the  object  can  be  estimated  directly  from  the  projected 
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(B)  (C) 


Figure  2.  Three  views  of  the  support  of  a  triangular  curved  surface.  (A)  x-z  view  of  the  sur¬ 
face.  From  this  view  we  see  that  the  height  distribution  is  the  same  for  every  value 
of  y  and  is  symmetric  about  the  x=8  plane;  (B)  x-y  view  of  the  surface;  (C)  y-z  view 
of  the  surface. 
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(C) 


(D) 


Figure  3.  Reconstruction  of  triangular  curved  surface.  (A)  x-y  support  bound.  The  bound  was 
obtained  by  intersecting  the  estimated  x-z  support  bound  [same  as  shown  in  Figure 
2(A)]  with  the  x-y  autocorrelation  support;  (B)  Magnitude  of  the  true  reflectivity; 
(C)  Magnitude  of  the  initial  guess  on  the  support  of  (A);  (D)  Magnitude  of  recon¬ 
structed  reflectivity. 
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autocorrelation  in  practical  scenarios  involving  measurement  noise  and  2-D  height  variation  (instead 
of  the  1-D  variation  of  this  example),  this  example  demonstrates  that  a  filled  surface  reflectivity 
can  be  estimated  with  a  loose  angle- angle  support. 

6.  CONCLUSIONS 

We  have  reported  the  first  demonstration  of  three-dimensional  phase  retrieval  with  an  opacity 
constraint.  These  results  make  feasible  the  possibility  of  an  entirely  new  three-dimensional  imaging 
modality.  Further  experimentation  is  required  to  determine  how  loose  the  height  constraint  can  be 
made  and  the  importance  of  the  role  of  the  initial  height  estimate.  We  believe  that  the  theory  of 
two-dimensional  support  bounds  [7]  can  be  generalized  to  the  three-dimensional  case  to  give  good 
support  bounds.  We  also  conjecture  that  opacity  will  play  a  role  to  help  tighten  range-dimension 
support  bounds.  Historically  modes  of  stagnation  have  been  studied  in  the  two-dimensional  phase- 
retrieval  problem  and  appropriate  algorithms  to  avoid  such  modes  have  been  developed.  A  similar 
process  needs  to  be  undertaken  in  the  three-dimensional  case.  Finally,  the  sensitivity  of  three- 
dimensional  phase  retrieval  to  noise  needs  to  be  investigated. 
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Abstract 

We  discuss  a  3-D  imaging  modality  called  Pulsed  Heterodyne-array  Imaging  (PHI).  The 
relationship  between  PHI  and  stepped-frequency  methods  for  3-D  coherent  image  formation  is 
derived.  For  both  cases,  we  consider  flood-illumination  of  the  object  and  detection  with  a  2-D 
array  of  coherent  receivers  located  in  the  pupil  plane.  It  is  shown  that  PHI  can  recover  the  same 
coherent  3-D  image  as  with  a  stepped-frequency  method  such  as  holographic  laser  radar  [1]. 

1.0  Introduction 

In  previous  research,  we  have  demonstrated  the  ability  to  form  3-D  coherent  images  by 
flood-illuminating  an  object  with  a  tunable  laser  and  detecting  the  scattered  light  with  a  2-D  array 
of  coherent  (field-recording)  receivers.  As  the  source  frequency  is  tuned,  the  optical  field  is 
recorded  and  a  3-D  image  is  formed  by  3-D  Fourier  transformation  of  the  data.  We  refer  to  this 
method  as  Holographic  Laser  Radar  (HLR). 

The  basis  for  HLR  imaging  is  the  3-D  Fourier  Transform  relationship  between  the  optical 
field,  U,  and  the  object’s  3-D  complex  reflectivity  function,  a(r).  This  relationship  can  be  written 
as  [2-4] 

U(rd,  co)  =  a  exp{ikR)  J  a(r)  exp  [-zoo  (p  •  r)  ]  d3r  ,  (1) 

where  rd  is  the  unit  vector  corresponding  to  the  location  of  a  detector  element  with  respect  to  the 
object  domain  origin  ( r  -  0),  co  is  the  source  angular  frequency,  a  is  an  inessential  complex  con¬ 
stant,  R  is  the  distance  from  the  source  to  the  object  domain  origin  to  the  detector,  and  p  =  (rs  + 
rd)/c,  where  c  is  the  speed  of  light  and  rs  is  the  unit  vector  corresponding  to  the  source  location 
with  respect  to  r  =  0.  For  HLR,  data  is  recorded  for  a  series  of  rd  values  as  the  source  frequency  is 
tuned.  For  this  reason  Eq.  (1)  is  written  as  a  function  of  these  two  variables,  rd  and  ca 

By  using  Eq.  (1),  we  have  been  able  to  demonstrate  fine-resolution  3-D  coherent  imaging 
[5],  For  these  experiments  we  used  a  tunable  dye  laser  and  a  256  X  256  coherent  detector  array.  A 
3-D  DFT  was  used  for  image  formation  and  range  resolution  on  the  order  of  10pm  was  reported. 
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2.0  Pulsed  Heterodvne-Arrav  Imagine 


Now  consider  3-D  image  formation  with  pulsed  illumination.  This  imaging  modality  has 
been  discussed  previously  [6];  however,  in  this  paper  we  present  the  underlying  theory  for  image 
formation  and  analyze  it  in  the  context  of  stepped-frequency  imaging.  To  begin  this  analysis  con¬ 
sider  the  output  of  a  mode-locked  laser.  A  mode-locked  laser  coherently  phases  longitudinal  cav¬ 
ity  modes  of  the  laser  to  produce  a  series  of  short  pulses.  In  frequency  space  the  output  can  be 
written  as 


A  (co)  =  [comb  (aco)  •  rect  (hco)]  <8>  8  (co-co0) ,  (2) 

where  0  represents  a  convolution  and  8  represents  the  Dirac  delta  function.  To  simplify  the 
expression  we  have  assumed  that  the  modes  all  have  the  same  phase,  however,  this  treatment 
extends  readily  to  the  more  general  case  of  linear  mode  phasing.  An  illustration  of  the  laser  output 
in  frequency  space  is  shown  in  Fig.  1;  the  pulse  consists  of  several  laser  cavity  modes  that  are 
equally  spaced  in  frequency. 
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Figure  1.  Frequency  Representation  of  a  Mode-Locked  Laser  Pulse. 


The  temporal  form  of  this  laser  output,  found  by  taking  the  Fourier  transform  of  Eq.  (1),  is 

*  (,)  =  ® ”C(D] ' exp ■  (3) 

which  is  a  time  series  of  sine  pulses  of  light  at  base  angular  frequency  co0- 

Now  consider  the  return  obtained  by  illuminating  an  object  with  such  a  mode-locked 
pulse.  It  follows  that  each  frequency  is  multiplied  by  a  complex  constant  that  is  equivalent  to  the 
detected  optical  field  that  results  when  the  object  is  illuminated  with  the  particular  laser  frequency 
at  unit  amplitude.  Using  Eq.  (1),  we  can  write  the  return  pulse  as  a  superposition  of  the  output  fre¬ 
quencies  with  each  multiplied  by  a  complex  constant  corresponding  to  the  return  field 

Aref(rd,  co)  =  U(rd,  co)  {  [comb  (aco)  •  rect  [boa)]  ®  5(co-  co0) }  ,  (4) 
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where  U(rd, co)  is  the  optical  field  from  the  object  for  the  particular  detector  location  and  source 
frequency  given  in  Eq.  (1).  As  shown  by  this  equation  and  illustrated  in  Fig.  2,  frequency  analysis 
of  the  return  pulse  yields  a  sampled  version  of  the  objects  complex  reflectivity. 
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Figure  2.  Frequency  Representation  of  a  Returned  Mode-Locked  Laser  Pulse. 


The  detection  scheme  for  PHI  is  to  detect  the  return  given  in  Eq.  (4)  using  a  heterodyne  receiver. 
The  resulting  output  is  the  optical  field  as  a  function  of  time.  This  temporal  signal  is  found  by  tak¬ 
ing  the  Fourier  transform  of  Eq.  (4)  with  respect  to  the  time  variable  which  gives 


(5) 


As  shown  in  Eq.  (4),  PHI  provides  a  convenient  means  for  sampling  the  object  field  at  a  series  of 
equally  spaced  samples.  Also,  2-D  Fourier  transformation  of  the  received  signal  for  each  time- 
bin  with  respect  to  the  variable  rd  will  yield  an  estimate  of  the  3-D  object. 


A  block  diagram  for  a  PHI  system  is  shown  in  Fig.  3.  Essential  elements  include  a  pulsed 
laser  source  and  a  2-D  array  of  pupil-plane  heterodyne  receiver  elements.  3-D  image  recovery  is 
achieved  by  2-D  Fourier  transformation  of  the  received  data  in  the  dimensions  of  the  receiver 
array  for  each  time  bin.  Note  that  only  the  return  from  a  single  pulse  is  needed  for  image  recovery. 

Further  insight  into  PHI  can  be  gained  by  considering  the  temporal  Fourier  transform  of 
the  received  field.  This  is  equivalent  to  the  received  field  when  the  object  is  illiminated  by  an 
impulse.  From  Eq.  (1)  we  have 

U(rd,t)  =  a  Jo(r)5(r-^  +  (^^)-r)d3r  .  (6) 

For  the  monostatic  case  with  rs  and  rd  pointing  in  the  z  direction,  we  have 

U(rd,t )  =  a  J  o(x,y,  z)5^r-^  +  2-^dxdydz  ,  (7) 

where  we  have  switched  to  Cartesian  coordinates.  Further  evaluation  of  this  equation  yields 
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Figure  3.  Block  diagram  for  Pulsed  Heterodyne-array  Imaging  (PHI). 


U(rd,t )  =  a  J  <5p(z)h[t-^  +  2^dz  ,  (8) 

where  op(z)  is  the  object’s  complex  reflectivity  projected  onto  the  z-axis.  Evaluation  of  this  inte¬ 
gral  yields 


&(r„«)  =*<*,(§-!')  (9) 

where  a’  is  an  innessential  constant.  Now  consider  the  received  field  in  the  time  domain.  From 
Eqs.  (5)  and  (9)  we  have  in  the  case  of  monostatic  operation  along  the  z-axis 


xret 
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exp  ( — i  cop?)  .  (10) 


This  result  shows,  in  agreement  with  intuition,  that  the  return  consists  of  a  time  series  of  sine 
pulses  convolved  with  the  projected  cross-section  of  the  object’s  complex  reflectivity. 


3.0  Summary 


In  summary,  we  have  shown  that  PHI  provides  a  way  to  generate  coherent  3-D  images  of 
targets.  It  uses  a  pulsed  laser  with  pulse-length  equal  to  the  desired  range  resolution.  The  use  of  a 
pulsed  laser  eliminates  the  need  for  performing  a  Fourier  transform  in  the  source  frequency 
dimension  as  required  with  HLR.  It  is  anticipated  that  PHI  will  be  useful  for  3-D  imaging  distant 
objects  where  pulsed  illumination  and  high-speed  detection  are  advantageous. 
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